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Abstract 



We study the decay of the inflaton in a general Z2 x Z2 symmetric two 
scalar theory. Since the dynamics of the system is dominated by states with 
large occupation numbers which admit a semiclassical description, the decay 
can be studied by solving the classical equations of motion on the lattice. Of 
particular interest is the case when the cross-coupling between the inflaton 
and the second scalar field is negative, which is naturally allowed in many 
realistic models. While the inflaton decays via parametric resonance in the 
positive coupling case we find that for negative coupling there is a new mech- 
anism of particle production which we call negative coupling instability. Due 
to this new mechanism the variances of the fields grow significantly larger 
before the production is shut off by the backreaction of the created particles, 
which could have important consequences for symmetry restoration by non- 
thermal phase transitions. We also find that heavy particles are produced 
much more efficiently with negative coupling, which is of prime importance 
for GUT baryogenesis. Using a simple toy model for baryogenesis and the 
results of our lattice simulations we show that for natural values of the cross- 
coupling enough lO^^GeV bosons are created to produce a baryon to entropy 
ratio consistent with observation. This is to be contrasted with the situation 
for positive coupling, where the value of the cross-coupling required to pro- 
duce such massive particles is unnaturally large. In addition to our numerical 
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results we obtain analytical estimates for the maximum variances of the fields 
in an expanding universe for all cases of interest: massive and massless in- 
flaton, positive and negative cross-coupling, with and without significant self 
interactions for the second field. 



Typeset using REVT^ 



I. INTRODUCTION 



It has recently been realized |]T|,H that an oscillating scalar field can transfer its energy 
efficiently to other bosonic degrees of freedom via a process analogous to the classical 
phenomenon of parametric resonance. The most important application of this mechanism 
is to reheating after infiation. When the infiationary stage terminates the universe is es- 
sentially devoid of matter and the infiaton field oscillates with large amplitude about the 
minimum of its potential. Its coupling to other bosonic fields can then lead to explosive par- 
ticle production via the non-perturbative resonance mechanism. This is sometimes called 
preheating in the literature because the produced particles are far from thermal equilibrium. 
The simplest model for studying this phenomenon is a theory with two real scalar fields and 
effective potential 



V{^, X) = + + + ^X' + f (1) 



Such models exhibit chaotic infiation when one of the fields, called the inflaton, acquires 
a large expectation value (compared to the Planck mass) 0]. The energy density of the 
universe then becomes dominated by the potential energy of the infiaton and the universe 
starts infiating. This continues until the amplitude of the infiaton expectation value becomes 
of order the Planck mass. At this point the infiationary stage terminates and the infiaton 
starts oscillating. 

An intuitive picture of the parametric resonance mechanism emerges by assuming that 
the infiation field oscillates sinusoidally about its minimum 0. With t) $o cos{uj'^t) 
and neglecting the A^^ term, the mode equations for the % field can be written as the Mathieu 
equation 

^ + [Aoik) - 2go cos(2^)]xfc = , (2) 



where z = u'^t and uj^{ky = k'^+m^ is the frequency squared of Xk- It is well known that the 
Mathieu equation possesses unstable solutions for which the modes grow as Xk oc exp[/ifcCo'°t]. 
This corresponds to exponentially growing occupation numbers n^{k) oc exp[2fikuj^t] and is 
interpreted as particle production. The most important features of the solutions to Eq. (0) 
are readily understood from the stability chart depicted in figures 1(a) and (b). The fi = 
curves divide the chart into stable (dark) and unstable (light) regions. Some curves of 
constant positive /i are also shown. For clarity figure 1(a) contains only a limited range of 
Aq and go, although the physically interesting values are |go| ^ 1 BHH- The go ^ 1 regime 
is called "broad resonance" in the literature. Figure 1(b) shows a piece of the instability 
chart for large |go| around the line Aq = 2go. Notice that, if one makes the usual assumption 
g > 0, then for > the allowed Aq lie above the line Aq = 2go (c/. Eq. (|)). As a 
consequence the regions of large /i are excluded, and the maximum value fi reaches is 1, 
attained for Aq ^ 2go and go ^ 1. Above Aq = 2go, /i is a rapidly decreasing function of 

1/2 

Aq. Since the separation between neighboring instability bands for large go is of order Qq , 
this implies that the dominant production of x particles occurs for k ~ /cj-cs — %^'^^d>- 
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is important to note that the inflaton also decays into its own fluctuations via parametric 
resonance 0. However, the value of /z corresponding to this decay channel is typically a 
few times smaller than the value of /i for decay into x fluctuations as long as go is large. 
Hence, unless the decay into x particles is somehow shut off, the decay into particles is 
sub dominant. 

While the above picture is simple and intuitive, in reality the situation is much more 
complex. First of all, the equations of motion (EOM), and hence the dynamics, are strongly 
modified in an expanding universe. In particular, the amplitude of the inflaton oscillations 
decreases even in the absence of particle production. Secondly, as the inflaton decays, 
the backreaction of created particles alters the parameters in the Mathieu equation. The 
amplitude decreases as energy is drained away and the masses get contributions of the form 
Sml = g{{5xf) + 3A^((50)2), 5ml = g{{5ct>f) + ^\{{5xf). with ((50)^) = (0^) _ (0)2 and 
((^x)^) = (x^) ~ (x)^ denoting the variances of the fields. One can still gain insight into the 
decay by replacing the initial values Aq and go in Eq. (^) by time dependent parameters 

+ (4) 



where 



X eff - '"x ^ y \v"v^; / "1" '^'^x\ 
ujl^^ = ujl + g{{5xf) + ^Mm'). (5) 

Here $ is the slowly varying amplitude of the oscillating inflaton expectation value 0o(t) = 
{(j){x,t)), and is the frequency of 0o(t) in the absence of backreaction. We distinguish 
LU^ from its initial value a;^ because due to the 0^ term in Eq. (P the frequency is time 
dependent even without decay in an expanding universe. The reason it is useful to think 
about the decay in terms of the parameters defined in Eq. @ is that the basic features of 
the resonance mechanism are extremely robust. The amplification of the modes depends 
only on the fact that one is dealing with oscillators with a time dependent frequency, and 
in what follows we will only make use of the most prominent features of the stability chart 
when obtaining analytical estimates. For large g these are independent of the details of the 
potential and the exact time dependence of the parameters. The description in terms of 
an equation of the form of Eq. (0) finally breaks down completely when the resonant mode 
amplitudes grow large. For then the nonlinear terms in the EOM become important and 
scatterings become fast. The result is a slowly evolving "scattering regime" characterized 
by smooth power spectra ||TT|-[13|. The resonant decay of the inflaton has been treated in 



various approximations by many authors |[I|,0,^-TD]. The full nonlinear problem has been 
studied in detail using lattice simulations. The idea is that one can evolve the system using 
the classical equations of motion because the dynamics is dominated by states with large 
occupation numbers which admit a semiclassical description. This approach was pioneered 
in where it was applied to the decay of a massless field into its own fluctuations. The 



decay into a second massless field with and without self-coupling was treated in |12], and 
the decay of a massive inflaton was discussed in ||T3 |. 



Our main interest in this paper is to investigate what happens when the coupling g in 
Eq. (P is negative, which has not been previously studied. The point is that parameters 
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like g that couple two or more scalar fields appear in almost all extensions of the standard 
model, and that there is generally no reason whatsoever that these parameters should be 
positive. Of course the effective potential should be bounded from below in order to ensure 
the existence of a sensible vacuum. For the potential given in Eq. (0) this requires 

r^^>l. (6) 

We will assume that this stability bound is satisfied. 

Naively, as one can see from figure 1(a), the negative g case is dramatically different 
from the positive g case. The resonant momenta are now limited to the region above the 



line Aq = — 2|go|, instead of above the line Aq = 2\qo\. When 2|go| > |^o| ^ 1 [Hll 



cosh 27r/i = cosh 3^ 



2tt 





X COS 3? 







dz{Ao - 2gocos(2z))2 
dz{AQ — 2qo cos(2z))5 



(7) 



Along Aq = — 2|go| this evaluates to /i = (4/7r)|go|"'^^^ (see also [T^). For large |go| this is 
easily orders of magnitude larger than the corresponding positive g value, which is always 
^ 1. Since the decay time is proportional to one might hence expect that particle 
production and consequently inflaton decay is much faster for negative g. 

This conclusion, however, is premature since we have implicitly assumed that initially 
the expectation value xo = {x{^^ t)) = 0- To see whether this is the case one has to examine 
the EOM for the fields during inflation. Note that the large expectation value (j)o{t) of the 
inflaton induces a 0o-dependent minimum in the potential for the second field: 

Xl = \'^ for ml + g^l<0 (g) 
1 otherwise 

It is not hard to see analytically that if xo ~ Xo initially, then it remains so during inflation, 
as long as \g\/X(i) ^ 1/2 and is not too large jl^. Our numerical simulations show that 
even if xo is not close to Xo initially, it will become so by the end of inflation, provided the 
above conditions are satisfied. Since the physically important broad resonance case with 

^ is our main concern, this means that the correct initial condition for the oscillatory 
regime is given by Eq. (|^). 

To illustrate the effect of Xo 0; consider the linearized mode equation for 6xk, the 
Fourier transform of x ~ Xo- Neglecting the expansion of the Universe we obtain 



cf5xk 
dH 



+ [e + ml + 3A^Xo + 9<i>i) SXk + 2g^oXoHk = (9) 



where 0o is the oscillating expectation value of the inflaton field. Note that if Xo ~ Xo 
even during the oscillatory regime, then the term in parentheses of (^) is positive (c/. (^). 
This would mean that the evolution of 6xk resembles that of a positive g case. Neglecting 
the last term for the moment, we see that in the simple case when = and xo = Xo, 
Eq. (I) reduces to Eq. (D with Aq = (fcV(tu^)2) + 2qo and go = 2|^|$g/(2w^)2. In the case 
7^ there is no such simple correspondence, but the effective q is still positive. Note 
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that the effect of the last term in Eq. (|^) is that it produces resonant growth of 5xk when 
the amphtude of dipk is large and Xo 7^ 0. Since most of the time 5xk ^ ^'Pk for the resonant 
modes, the effect of this term will be small. The analogous term in the mode equation for 
speeds up its growth when 5xk is large. 

We have just seen that if xo follows its instantaneous equilibrium value xo the fact that 
g is negative should not make much difference. The situation changes drastically if xo is not 
"enslaved" by 0o in this way: there is then a portion of the period for which the term in 
parentheses of is negative. During this time the 5xk are unstable and are expected to 
grow exponentially with large /x, just as in the naive negative g picture discussed initially 
(c/. also [|15|). The physics hence depends crucially on the detailed dynamics of the field 
expectation values. This dynamics is very complicated due to the interactions with the 
nonzero modes, especially once their amplitudes become large due to the resonance. To 
fully investigate the negative g case one must hence resort to lattice simulations, and this is 
the subject we turn to next. 

In section |I| we will discuss massless fields, followed by a detailed study of massive fields 
in section |ITT[ The relevance of our results for grand unified theory (GUT) baryogenesis are 
discussed in section |^ Finally, our conclusions are presented in section 0. 

Before diving into the analysis we should point out that most of our numerical calcu- 
lations and analytical estimates will be concerned with obtaining the maximum variances 
attained by the fields during preheating. There are two reasons for focusing on these quan- 
tities: first, the maximum variances determine the symmetry restoration "power" relevant 



to the occurrence of so-called nonthermal phase transitions |ll8|JT9|| . Second, for very massive 



particles such as the bosons needed for GUT baryogenesis the variance is directly propor- 
tional to the energy and particle number densities. This will be discussed in more detail in 
sections fT^ and |V|. 



II. MASSLESS FIELDS 
A. Numerical Results 

We begin our investigation by studying massless fields, which is convenient for two rea- 



sons. First, the positive g massless case has been studied extensively on the lattice in [[12 
allowing for a detailed comparison of the positive and negative g cases. Second, massless 
fields (m^ = = 0) are the simplest to treat numerically because in this case the EOM 
in a Friedman- Robert son- Walker (FRW) universe are conformally equivalent to those in 
Minkowski space. By this we mean the following: in terms of the variables r = / dt/a{t), 
= 0a(r)/a(O), and x = X'2(''")/a-(0), where a(r) is the scale factor, the expanding universe 



EOM reduce to those of two interacting scalar fields in static space-time ||Tl|JT^ . This is 
true provided the energy density behaves as radiation dominated, which is indeed the case 
for a massless inflaton in the oscillatory regime 0. Since it is less time consuming to do 
simulations in Minkowski space our numerics for the massless fields were done in terms of 
the barred variables defined above, and our lattice time was r. The conversion of lattice 
quantities to physical quantities will be discussed in detail below. 



The computational techniques involved are described in detail in In short we evolve 
the classical EOM on a three dimensional lattice starting from certain random initial con- 
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dition for the inhomogeneous modes |T^. The initial conditions for the zero modes (field 
expectation values) are given by the values at the breakdown of the slow roll inflationary 
stage: $0 = $(t = 0) ~ Mp, where Mp = (STrG)-^/^ ^ 2.4 x lO^^GeV is the reduced Planck 
mass, and Xo = Xo as in Eq. (||). The results presented in this section were obtained on 128^ 
lattices. 

We begin by comparing the dynamics of two theories with opposite signs of g but oth- 
erwise identical parameters. The chosen values are = 10~^^, A^^, = 10~^, and \g\ = 10~^°. 

1/2 

For a massless infiaton the zero mode oscillates with frequency uOip = cA^ $, where 
c ~ 0.85. Hence these couplings correspond to go — 35 and r = 10 (c/. Eqs. (^), (^). 
The evolution of the expectation values 0o and xo for the negative and positive g cases are 
shown in figures 2 and 3, respectively. The corresponding variances are shown in figures 4 
and 5. Since the energy density in the infiaton zeromode is proportional to we see from 
figure 2 that in the negative g case 40% of the initial energy is in fluctuations by the time 
r ^ 3 X 10^$o ^- Notice also that Xo follows 0o faithfully {i.e. Xo ~ Xo) until r pa 8 x 10'^$o \ 
After this the two expectation values are no longer in phase, even though the dynamics of xo 
is clearly still strongly influenced by (po. Figure 3 shows that it takes until r ^ 4.5 x 10*$q ^ 
for 40% of the energy to decay in the positive g case, and in contrast to the situation for 
g < here almost all of the decayed energy goes into (f> fluctuations. After r ~ 4.5 x 10®$o ^ 
the decay rate becomes very small, a fact that is also reflected in the variances (figure 5), 
which are both changing very slowly beyond this point. Similarly we see in figure 4 that for 
negative g the variances level off at r ~ 9 X 10^$o \ although the variance is still growing 
slowly until r ~ 2.2 x 10^$o ^. We will explain all of these features, and give estimates for 
the magnitudes of the relevant quantities, in section 11 hi . 

Before doing so let us briefly discuss the conversion of lattice quantities to physical 
quantities. As discussed at the beginning of this section, the lattice fields are rescaled 
by a factor a('r), where r is the conformal time (for convenience we define t = r = 
to correspond to the beginning of our simulation (i.e. the end of slow roll inflation) and 
choose a(0) = 1). This means that one must simply divide lattice quantities such as energy 
densities and variances by the appropriate powers of a(r) to get the physical values. For 
example, (Sx^) = {^X^) / di^)'^ ■ In the oscillatory regime of a massless infiaton the scale 
factor behaves as in a radiation dominated universe, so we can write a(r) = {I + Hqt) where 
Hq = (p(0)/3)^/^ Mp ^ is the Hubble constant at the end of inflation. It is easy to show that 
at the end of slow roll $0 ~ 1.8Mp and the kinetic energy density is about half the potential 
energy density |2^. Using the fact that ppe = A(^$q/4 foYg>0 and ppe = (1 — ^/^)^<f>^o/^ 



hi g <0 we thus obtain Hq ^ V^^Mp for 51 > and Hq ^ {X^{1 - l/r)y^'^Mp for g <0. 

With these estimates, and recalling that our lattice time is r, it is simple to convert 
from lattice to physical quantities. This was done in obtaining figure 6(b) from figure 6(a): 
figure 6(a) shows the maximum lattice variances for negative g as a function of go for two 
values of r. The corresponding physical variances are shown in figure 6(b). Since the 
variances are fiuctuating significantly (see, for example, figure 4), both peak and valley 
values are given. These curves, too, will be explained in section [11 B[ Here we simply note 
that when the lattice variances behave as shown in figures 4 and 5 it is easy to obtain the 
maximum physical variances: since we must divide by a(T)^ to obtain the latter from the 
former, the maximum occurs at the time when the slowly varying state sets in. After that 
the slow increase cannot keep up with the decrease due to the expansion. 
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B. Discussion and Analytical Estimates 



The most striking feature of figures 2-5 is the fact that after an initial period of rapid 
exponential decay the system reaches a slowly evolving state during which the rate of energy 
transfer from the inflaton zeromode to the fluctuations is quite small. This effect has been 
studied in detail in |T^, where it is shown that for large go or r backreaction and scattering 
limits the maximum amplitudes of the variances and essentially shuts down the resonance. 
Throughout this paper we will refer to the slowly evolving state as the scattering regime 



21| . The variances of the physical fields at the time the scattering regime begins are of 
central interest since they are the maximum variances of the fields reached during preheating. 
The reason that the decay of the inflaton shuts off for large go or r while most of the 
energy is still in the oscillating zeromode is the backreaction of the created particles on 
the EOM. In particular the effective masses of the fields get contributions of the form 
5ml = dii^xy) + 3A<^((50)^), 5m^ = g{{5(j)f) + 3A^((5x)^)- As was explained in the 
introduction, the instability index /i^ is a rapidly decreasing function of A above the line 



A = 2\q\ (see also figure 1). In fact the decay due to the second band above this line 
is already very slow. Now the backreaction terms add a contribution 6 A = Srn^/uj^ to A, 
and for large g the distance between the first two bands above A = 2|g| is about \q\^^'^ {cf. 
figure 1(b)). Hence when the shift in A becomes so large that 

A{k = 0) + 5A>2\q\ + \q\^ (10) 

the first band above A = 2|g| is rendered ineffective and the decay slows down dramatically 
From this condition one can obtain estimates of the variances at the beginning of the 



slowly varying scattering regime. As explained in |jT2[, the (lattice) variances are then kept 
roughly constant at these maximum values by a feedback mechanism. It is important to 
point out that scattering plays a crucial part in the shut off mechanism. For example, we 
will see below that when A^^ is small the resonance is shut off by the g{{5(j)Y) term. This 
term grows large due to scattering of resonant x particles off the zeromode, an effect that 
is completely lost in the Hartree type approximations often used to study preheating. The 
fact that the (f) variance is responsible for shutting off the resonance for small A^^ was first 



realized in [12 



Before deriving estimates for the variances based on Eq. (|I^) we wish to point out that 
there is always a substantial "intrinsic" uncertainty associated with such estimates. By 
this we mean that physical quantities such as maximum variances, energies, et cetera, are 
extremely sensitive functions of the parameters in the model. For example, we have seen 
above that for large g a shift of order |g|^/^ is enough to move the initial position of the first 
resonant momentum above A = 2|g| from fc^ to /c^ ~ Since fik is a rapidly 

decreasing function of k p2| such a shift can have a substantial effect on the early evolution 
of the system. The upshot is that for a small change 5g/g = Sq/q ~ |g|^^^^ ^ 1, the time 
at which the scattering regime sets in changes by a factor of a few. As a consequence the 
maximum physical variances reached vary by about an order of magnitude. In this context 
we point out that the maximum lattice variances reached are much less sensitive to small 
changes in the parameters. The reason is simply that the lattice variances do not decrease 
due to the expansion, so the time it takes to reach the scattering regime is immaterial. This 
observation explains why the curves in figure 6(a) are relatively smooth compared to those 
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in figure 6(b). The sensitivity to small changes in parameter values has also been observed 
ini. 

We will now estimate the maximum variances in the massless case based on Eq. (|10|). 
Positive and negative coupling - Cases I and II - will be treated separately. As explained 
above, the maximum physical variances are reached at the beginning of the slowly evolving 
scattering regime. Since the variances fluctuate we distinguish between peak values and 
valley values. Estimates obtained from Eq. (|TUp correspond to peak values since Eq. ( pUj ) is 
the condition for shutting off the resonance. This condition must be met during those parts 
of the period of the infiaton zero mode when the resonance is active, but in between the 
variances can drop to smaller (valley) values. The formulae below will be given in terms of 
the quantity ^s, which is the amplitude of oscillations at the time the maximum variances 
are reached p3. We will discuss how to estimate $s after deriving the expressions for the 



variances. Note that Eq. (0) involves the time dependent parameter q rather than the 
initial value go- Recalling from section [II A| that in the massless case u;^ = c\r^, with 



c ~ 0.85, a glance at Eqs. (^ and (^ reveals that both the numerator and the denominator 
of q scale as a(t)~^. Hence, while q may change significantly due to the decay of the infiaton, 
it does not change due to the expansion in the massless case. 
Case I. For g > 0, Eq. (|10D reduces to 



6A = — ~ g2 , (11) 

^<t> eS 

with "^^gfj and uj'^es given by Eq. (||). When the self-coupling is small, in the sense that 
the 3Xy^{{5xf) term in fn?^ may be neglected compared to g{[S(f)'^)), Eq. ( |1 1|) yields an 
estimate for the field variance: 

((50)')pcak ^ ^^"^$s^0 cff = ^ g"^$s • (12) 

This can be turned into an estimate for ((5x)^) by assuming approximate equipartition 
of energy between the fluctuations of the two fields. In ||12| we have argued that this is 
reasonable because the transfer of energy between the two fields is efficient in the scattering 
regime. We thus obtain 

{i6x)\e.,^f^ = ^q-'^l (13) 
Our original assumption that was small can be translated into the constraint 

{m')<^,q-'^r-'^!. (14) 

Hence Eqs. (|T^ and (|TB|) are valid provided that r < (3c^g^/^)~^. This condition defines 
what is meant by "small A-^". Note that it translates roug hly into A^ < g^^^Xl^^ None of 
the runs presented in this paper satisfy this condition, and the above estimates are included 
here only for completeness. The "small A^^" case was treated in detail in |T^. 

For "large A^^", i.e. rq^^"^ > (3c^)~^, Eq. (|ll|) immediately yields an estimate for the x 
field variance: 
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((<5X) )pcak ^ (15) 
^ g-ir-i$2. (16) 



48c' 

This is obtained using tu^gg ~ c^;?, which is true because here r > (12c2gi/2)-i, The run 
shown in figures 3 and 5 belongs in the large A^^, positive g category. To see how well 
Eq. ([I6| ) works in this case note that we can replace {SxY with the lattice quantity on 
the left hand side if we do the same with $s on the right hand side. From figure 3, $s ~ 
so Eq. (|TB|) predicts (((5x)^)peak = 1-4 x 10^^$q. From figure 5 we see that this is in good 
agreement with the first peak ((5x)^)peak — 

1 X 10-^$^ at r ^ 0.65 x 10^$o\ To test the 
qo dependence in Eq. (|TB|) we ran our code for several values of go ranging from 3.5 to 1000. 
A power law fit to the data gives ((5x)^)poak oc qo^'"^^, in reasonable agreement with the 
prediction. 

Figure 5 also illustrates nicely that the infiaton decays into its own fiuctuations via 
parametric resonance, as discussed in the introduction. Once the x resonance is shut off 
this usually subdominant process becomes important and the variance continues to grow. 
After r = 2 X 10^$o ^ the growth is clearly exponential, and during this stage the inflaton 
decays essentially as if it were not coupled to the x field at all. The decay finally stops when 
the field variance reaches its scattering value ((50)^)peak ~ 10^^^'^{Tf) ~ 10"^M| p5 



This is discussed in more detail in |]T2[, as is the amount of energy in the x fiuctuations at 



the beginning of scattering regime ||26 . 



Case II. When g < 0, Eq. (0) reduces to 



2 

5A=^~4|g|. (17) 



eff 



Since we must now allow for the possibility that Xo has an appreciable amplitude, Eq. 
for and u"^ should really be modified to 

^eff = ^7((50)') + 3A,(((<5x)') + x9 , 

u^l+9{{isxr)+xi)+^\,mr). m 



In this analysis, however, we set xo to zero. The reason is simply that we are interested in 
obtaining an estimate of the peak variance. Once ((^x)^) becomes appreciable xo longer 
follows Xo of Eq- ® [0]; and there are certainly times when xo is small (c/. figure 2). These 
times are precisely the times during which ((^x)^) reaches its maximum value, as can be 
seen from Eqs. (^) and (|TB|), as well as figure 4. 
With Xo = Eq. ([TtD simply says 

^leS^\9\^h (19) 

independent on u^. Assuming that {{Sx)^) dominates ^^cfr [EH obtain 

((^X)')pcak ^ JfL$^ = ^ \q\-'r-'<^l . (20) 
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This is a factor 4|g|^/^ times larger than the corresponding positive g value Eq. ([T6D . The 
difference between positive and negative g is nicely illustrated in figure 1. In order to shut 
off the fast exponential growth in the g > Q case, it is sufficient that the first instability 
band above A = 2\q\ is rendered ineffective, i.e. SA ~ g^/^. In the negative g case, on the 
other hand, the instability is shut off only when A is shifted all the way from A = —2\q\ to 
above the line A = 2\q\, i.e. 6 A ~ 4|g|. After ((5x)^)peak is reached, ((50)^) wants to grow 
to its one field scattering value ~ 0.1$^(rg*) ||2^, just as in the positive g case. However, due 
to the large x occupation numbers scatterings become fast before the (p variance can grow 
this large, and it is cut off at a somewhat smaller value. If one assumes that, because of the 



efficient scattering, equipartition is reasonably well satisfied, one obtains the estimate 29 



I 3 



Here we have taken the appropriate value of $ to be $s even though the (p variance really 
peaks somewhat later than the x variance. To get an order of magnitude estimate the 
difference should not matter as long as scatterings terminate the resonance soon after 
{{Sx)"^) reaches its maximum. While our simulations indicate that this formula holds fairly 
well for moderate values of q 100) we are at present unsure how well it extrapolates to 
very large q or r. The reason is that we cannot reliably capture both the and x resonances 
simultaneously if q is too large [0. We should point out that we expect Eq. (^) to fail for 



very large r (i.e. large A^). The reason is simply that in this case the x resonance shuts off 
for relatively small occupation numbers, so that (p ^ x scatterings should not affect the (p 
resonance very much, and hence (((50)^) should reach its one field value. 

Note that the variances in Eqs. (^) and (|2T|) do not significantly change uj^p at the 
time the scattering regime is reached since 5'((5x)^)peak5 3A0((50)^)pcak < is satisfied for 
\q\,r > 1. We can now compare the predictions of Eqs. ( pO|) and (0) with the numerical 
results shown in figure 4. Recall that we can turn the above estimates for physical variances 
into estimates for lattice variables by replacing the fields on both sides of the equation with 
the barred lattice fields. Then noting that $s ~ $o "we read from figure 4: ((^x)^)peak — 1 X 
10-3$2 and ((50)2)p,,k ^ 5 x lO'^^^. Eqs. @ and (|lD predict ((5x)')pcak = 0.33 x lO'^^^ 
and ((50)^)peak = 4 x 10~3$q, so that the agreement is quite good. 

We will now address the question how to estimate the quantity $s in terms of which our 
expressions for the variances, Eqs. ([T^), (]TB|), (|T3p, ([2D|), and (^), are given. Recall 

that $s is the amplitude of 0o at the time when the scattering regime is reached and the 
variances peak in the expanding universe. The first thing to note is that in the absence 
of expansion the amplitude of the zeromode does not decrease significantly by the time 
the scattering regime is reached. This can be observed in figures 2 and 3. The reason is 
simply that the energy density is proportional to so that even a substantial loss of energy 
corresponds to a rather small change in amplitude. To obtain an estimate for $s we can 
thus take the decrease to be due to the expansion alone. Then 

•^-s = ^ , (22) 

where Tg is the beginning of the scattering regime. As already pointed out at the beginning of 
this section, depends sensitively on the initial position of the resonance and is consequently 
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difficult to estimate accurately. Since the occupation numbers grow as n = no exp(2yUCo'°r) 
we find 

rs = ^ln(^), (23) 

where ngcatt is the occupation number of the resonant modes at the beginning of the scat- 
tering regime. Our initial conditions correspond to "rio"~ 1/2 |]12|, and nscatt is roughly 



for the small A;^, g > case, and 10 /X^ for the large A^, g > and g < cases |3T|. The 
difference arises because in the ffist case ((50)^) is responsible for shutting off the resonance 
while in the latter two cases the growth of ((^x)^) terminates the exponential regime (see 
also [|12|). This leaves us with the question what value to take for /i. From figure 1 we 
see that fi ~ 0.2 is reasonable for the positive g case, where we are limited to the region 
above A = 2\q\. This is also correct in the negative g case: it is not hard to show that 
Xo ~ Xo until {{Sx)"^) ~ (12c^r|g|^/^)~-'^$^ ( 10]), which is already larger than the positive g 
value given in Eq. (|T^). As discussed in the introduction, the negative g case is equivalent 
to positive g as long as xo ~ Xo, so during this initial period the positive g value for /i is 
appropriate. Once Xo stops following xo the variance rises rapidly to its negative g peak, 
in a time interval negligible compared to the already elapsed time. Hence we can estimate 



Tg using the positive g value for /i ||33|. For example, for the runs in figures 2-5 we obtain 
Ts ~ 1 X 10*$o ^5 reasonable agreement with the numerical results. Eqs. (UTj and (p3[) 



together with the formula a(r) = (I + Hqt) ^ (1 + [A(^(1 — l/r)]^/^Mpr) given in section [II A| , 
complete our estimate of $s. 

The final task of this section will be to explain the evolution of the variances in more 
detail. A closer look at figures 4 and 5 reveals a rather complex behavior of the variances. 
The first prominent feature is a slow modulation of the maximum variance. A simple 
explanation of this phenomenon can be gotten by assuming that the infrared (IR) modes 
oscillate in phase and can be modeled by one oscillator of a definite frequency uj"^ = k'^+uj'^ ^g. 
For -C tu^ eff we then have two weakly coupled oscillators, the IR mode and the zero mode, 
with the same natural frequency. These modes will transfer energy back and forth much 
like two pendulums coupled by a spring. As can be seen in figures 3 and 5, the maxima of 
the variance envelope correspond to minima of the zeromode envelope, and vice versa. Note 
that no such modulation occurs for the x fi^ld since no dominant x zeromode develops in 
the scattering regime. 

The second prominent feature of the variances is that in the scattering regime they fluctu- 
ate between "peak" and "valley" values with frequency 2Lj^eS- Understanding this behavior 
is potentially important for applications were precise values of the variances are required, 
such as nonthermal phase transitions and baryogenesis. A ffist estimate of the amplitude of 
these fluctuations has been given in Ref. [jl3[. Here we will present an alternative derivation. 



leading to somewhat different results. We will attempt to explain the origin of the variance 
fluctuations using a simple toy model which contains some generic features of the resonant 
growth, and also teaches us something about the slowly varying scattering regime. We will 
illustrate our model on the x fi^ld and ffist treat the case g > 0. After neglecting the non- 
linear (scattering) terms, the equation of motion for a field mode X = Xk can be written 
as 

X + uj\X = , (24) 
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where 



co'x = k' + gm + g{m') + 3A,((5x)') • (25) 

Since we are interested in IR modes that contribute significantly to the variance we can 
neglect the momentum dependence of the frequency. Using Eq. (|Tl|) we then find that in the 
scattering regime uj\ varies between cu^^x — — ^^^d cu^j^ ~ 3A^((5x)^) ~ 
g^/^C(j0. Dropping k"^ in Eq. ( P^D hence amounts to assuming that the variance is dominated 
by modes with k"^ < k"^^^ ~ q^^'^u^. Consider now the solution to Eq. ( ^i)) if ojx = for 
t < and uJx = for t > 0: 

X = Xq cos oj^t H sin uj^t , for t < 

UJQ 

^0 

X = Xq cos Uit H sin ujit , for t > . (26) 

UJi 

The important feature of Eq. (^) is that for uq > uJi the amplitude of the "kinetic term" 
Xo/a;i at t > is amplified by cuo/ix'i in comparison to the kinetic term at t < 0. This 
amplification captures the essence of the amplitude growth mechanism, and it may be used 
to explain the variance fiuctuations in the scattering regime, as follows. Consider a sequence 
of matchings as in Eq. (pBD , with Un/uJn+i = e^ Assuming that at each matching the 
solution of the previous time interval has a random phase pi] one obtains an amplification 
Xn/Xq = [(1 + e^'')/2]^/^. Using the fact that ojq/ujn = i^max/i^min = e'^*^ we find for e small 
that Xj^i/Xq ^ (u^inax/^^min)^^^- Siuce thesc are the modes that dominate the variance, this 
translates into ((5x)^)pcak/((5x)^)vaiiey ~ Xj^/^o ~ t^max/t^min- Plugging in the values for 
the frequencies obtained below Eq. (|25|) yields 



((5X) )pcak 



((5x)^)valley 



X max 



g3 ig>0). (27) 



This estimate agrees well with figure 5 after r ~ 4.5 x 10'^$q. Before that the oscillations 



are somewhat larger. Eq. ( pTf ) also agrees well with the large q runs presented in [0]. 

For negative g the fiuctuations in the variance are of different origin. For one thing, the 
peaks in iiow occur when 0o is maximum, in contrast to the positive g case just 

discussed. As already pointed out following Eq. (]T^), the largest peaks occur when xo ^ 
and 00 is large. The reason for this is simply that when xo gets out of phase with xo: the 
negative gcj)^ contribution to the x mass tends to destabilize the system. To compensate, 
the positive 3\^{{6xy)peak term must acquire a magnitude 3A^Xo5 since, as shown in the 
introduction, this is the value required to stabilize the mode equations. We have argued 
above that Xo falls out of phase with xq when ((^x)^) ~ (12c^r|gp/^)~^$^ [0. For example. 



in figure 2 this occurs at T ^ 9x 10^<l>o^ Fi gure 4 shows that the largest variance is reached 
a quarter of an oscillation later, when 0o reaches its maximum. During this first half-period 
after xo stops following xq there is genuine instability, causing the variance to grow. At later 
times the IR modes have sufficient amplitude to compensate the negative mass term without 
growing exponentially: they simply get pulled by 0o and keep the effective mass positive at 
all times. When (pQ passes through zero the variance returns to its previous "valley" value. 
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We can estimate the ratio of peak to valley variances by again appealing to our toy 
model of matching oscillatory solutions with different frequencies, as in Eq. (|^). There are, 
however, significant differences between the case at hand and the positive g case discussed 
earlier. For negative g the maximum frequency of an infrared x mode that contributes 
significantly to the variance occurs when 0o passes through zero and is approximately given 
by cu^ax — 3A^((5x)^)vaiiey This frcqucucy holds only for a small fraction of the oscillation: 
as soon as \g\(t>l > 3A;^((5x)^)vaiiey the x mode gets dragged along with the zeromode and 
its frequency changes to ciJmin — ^(i,es- Hence there is no continuous change in frequency as 
in the positive g case, and consequently only one matching of solutions a la Eq. (^). The 
single matching gives X^^,^/ X^i^ ^ ^^max/'^min, where we assumed that after the matching 
the kinetic term dominates. This is the case unless Xq/ujq ^ Xq in Eq. (^61) . Next we use 
((^X)^)pcak ~ Xo ~ I^I^V^x to obtain 

((^X)')pcak ^ /t^xmax\ | |i < Q) . (28) 

my) valley \ / 

This prediction can be compared to the numerical results presented in figure 6. As explained 
previously, the maximum lattice variances are much smoother functions of the parameters 
than the physical variances, and we will thus concentrate on figure 6(a). Note that even 
though we have presented |go| = 0.35 runs in figure 6, these runs just barely reach the scat- 
tering regime, especially for r = 2 (c/. 0]). After discarding the |g| = 0.35 runs we obtain 
as a best fit ((5x)^)peak ^ 0.1Q\q\-°-^^r-^-^^^^ and ((5x)^)vaiiey ^ 0.1\q\-^-^^r~°-^'^^^, and 
hence ((5x)^)peak/((^x)^)vaiiey — 1.6 1 g| This agrees reasonably well with Eqs. (pOD 
and (^). Note that the slopes in figure 6(a) are increasing with which may indicate that 
at |go| = 3.5 the variances still have not reached the full scattering regime values. This could 
be the explanation of the somewhat low slopes quoted above. Unfortunately we cannot test 
our estimates reliably for |g'o| > 350 since this would require enormous computing resources. 
The problem is essentially that as g and get large the particles scatter readily into high 
momentum states and one needs an enormous ultraviolet cutoff in order to accommodate 
the decayed energy. This must be combined with good IR resolution in order not to miss 
scatterings such as x{^res)<P{k = 0) — >■ x(^)0(^) which involve very IR momenta. The upshot 
is that we have found negative g simulations with |go| ~ 500 to be unreliable even on 128^ 
lattices (in the massless case). 



III. MASSIVE FIELDS 



A. Numerical Results 



For massive fields the expansion of the universe cannot be taken into account by a simple 



rescaling as was done for massless fields in section P A[ The equations of motion in a FRW 
universe are 



d'^x{x, t) ^ ^jj dxjx, t) ^ 



d(J)(x, t) 
ot 



dt 



a(t)2 
a(t)2 



+ 9X^ix,t) + A<^0^(f,t) 



+ 



\x,t) + \^x\x, t) 



0(x, t) = 
X(x, t) = 



(29) 
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where 



ait) V 3M| 



m-^-^ = J^, (30) 



and p is the total energy density. This set of partial differential equations was solved without 
approximation on the lattice. The initial conditions for the fields were chosen as follows: 



the (j) zero mode 0o was set equal to 2Mp. As explained in section [11 A| , this corresponds to 
a time slightly before the end of inflation, xo was set equal to Xo as given in Eq. (H), and 
the velocities <po and Xo were set equal to their slow roll values. The initial conditions for 
the inhomogeneous modes were chosen as described in |]T2[ and |]TB| . 



Here we will present a brief survey of our main numerical results for massive fields. 
In section p.11 B| below we will discuss the figures in detail and give analytical estimates for 



various quantities of interest. We should mention that throughout this paper we have chosen 
parameters that are "realistic" within the context of chaotic inflation. In particular we use 
$0 ~ Mp, Xfj, ~ 10~^^, and ~ lO^^GeV. The latter two values are determined by the 
observed anisotropies in the microwave background, as measured by the cosmic background 
explorer (COBE) |3|]. 



All of the numerical results presented are for g < 0. In figure 7(a) we show the maximum 
X variances produced in the expanding universe for a massless inflaton as a function of 
for various values of |go|- Figure 7(b) is a blow-up of the |go| = 350 curve. The important 
feature to observe is the extremely spiky nature of the curve, with huge jumps in production 
for small changes in the parameters. We have already discussed in section |II B| that the 
maximum variances are not smooth functions of the parameters, but the variation here is 
more extreme and of a different origin. It will be explained in section [111B| . Figures 8(a) 
and (b) show the variances and occupation numbers for a set of parameters for which 
there is very little growth, i.e. a point corresponding to one of the valleys in figure 7 ^7^. In 
figure 9 we show the same quantities for parameters were significant production does occur. 

Figure 10 shows the maximum peak and valley x variances produced in the expanding 
universe for a massive inflaton as a function of for |g'o| = 350. Again we observe extreme 
sensitivity to small changes in the parameters. In figures 11(a) and (b) we show the variances 
and X field zeromode as a function of time for three runs with just slightly different values of 
m~^. As can be seen in figure 11(a) the maximum variances reached differ by several orders 
of magnitude. In section |III B| we will explain how the growth of the variance depends on 



the dynamics of the x zeromode, which is illustrated in figure 11(b). 

An important point which will be discussed in detail in section [111B| is that a theory 



with negative coupling can produce massive particles for much smaller values of | go I- In 
fact we will see that for most of the parameter range in figures 7 and 10 there would 



be no particle production at all for g > 0. As will be explained in section the fact 
that massive particles are produced more readily with negative g is crucial in constructing 
GUT baryogenesis models with "natural" coupling constants. This will be illustrated by 
figures 12(a) and (b), where we show the maximum peak and valley variances reached for 
= 10^'^GeV as a function of \g\. For the range of couplings shown in the figures it is 



impossible to produce such heavy particles with g > 0. 



15 



B. Discussion and analytical estimates 



In this section we explain the main features of inflaton decay into massive particles, with 
particular emphasis on the resonant growth and its subsequent shut-off. We will also give 
analytical estimates for various quantities of interest. The physics of the shut-off mecha- 
nism is rather simple. For positive g parametric resonance produces particles with physical 
momenta of order kf^^^ ~ g^/^tu^, so that kf^^^ 3> when g ^ 1. When the particles 
mass exceeds the resonant momenta, the production shuts off. For comparison we note that 
perturbative decays are kinematically forbidden if the total mass of the decay products is 
greater than the energy of the decaying particles. As a consequence, parametric resonance 
can produce much heavier particles {m^ ~ q^^'^u^) than perturbative decays {m^ ~ u^). 
(When g -C 1 the resonant and perturbative scales coincide.) This observation could have 
important consequences in that parametric resonance provides a mechanism to produce the 
heavy gauge and Higgs bosons necessary for GUT baryogenesis. This will be discussed in 
some detail in section 0. Before presenting our results regarding the production of massive 
particles we wish to point out that from the outset one might expect the sign of the coupling 
g to play an important role in this context. The reason is that in the negative g case the 
maximum resonant momenta are of order kf^/^ ~ 2|g|^/^ti;0, which is by a factor 2|g|^/^ larger 
than in the positive g case. Therefore, one expects that with negative g the resonance is 
much more effective in producing massive particles. This expectation is confirmed below. 

For completeness we will study all cases of interest: massless and massive inflaton, posi- 
tive and negative cross-coupling g. We will then compare our estimates for positive g with 
the numerical results of 0, and for negative g with figures 7—12. The estimates we obtain in 
this section will be used in section ^ to estimate baryon production. In connection with this 
we remark that very massive particles are created only marginally relativistic and quickly 
become nonrelativistic due to redshifting. In this situation the variance is simply related to 
the number density of particles n and energy density p of the field: ~ m^((5x)^) and 

Based on Eq. (|l^) one can make a quantitative estimate of when the resonance shuts 
off. Indeed, when A becomes larger than about 2|g| + Iq]^^"^ the inflaton decay slows down 
dramatically, essentially because the instability exponent p decreases rapidly above the 
A = 2\q\ line (c/. figure 1). If the increase in A is mostly due to the (tree level) mass m^, we 
say that the resonance is shut off by the x mass. If, on the other hand, the growth of A can 
be attributed mainly to the backreaction effects (the contribution of the growing variances to 
the effective mass), the situation resembles the massless case in the sense that the variances 
reach their scattering regime values of section [II B| and a slowly varying state sets in. Before 
we begin studying the details of each of the cases mentioned above, we point out that the 
Universe expansion is crucial in shutting off the inflaton decay. This is so in essence because 
the inflaton amplitude decreases as the Universe expands, while the tree level mass stays 
constant. More concretely, the resonance shuts off when the mass becomes greater than the 
typical physical resonant momentum, which scales as k^^/^ oc a^^ for a massless inflaton 
and kP^/' oc a ' for a massive inflaton. This means that, as the Universe expands, the 
(physical) scale on which particles are produced redshifts and hence the relative importance 
of the X niass increases, leading eventually to the termination of particle production. Once 
the resonance is shut off, the only energy exchange mechanism that remains active is the 
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perturbative scatterings. We now turn to discuss the details of each of the cases mentioned 
above. 



1. Massless inflaton, g > 

For positive g and a massless inflaton the criterion for the x resonance shut-off is (c/. 
Eq. O) 

^eff^g^eff-glK)'(7)'' (31) 

where m^gg = + + 3A^((5x)^)- Here we have exphcitly written the dependence 

on the scale parameter a, used the fact that $ oc for a massless inflaton, and assumed 
ujfPeS — ^<j>- For massless fields the maximum variances reached were given in Eqs. (|T2p 
and (|13D for "small" (rgV2 < (Sc^)-!) and in Eqs. (|15D and (|16|) for "large" A^. It is 
clear that these estimates still hold for a massive x fi^ld if, by the conformal time Tg at 
which the scattering regime is reached, < g^/^(ci;°)^(ao/a(Ts))^ is satisfied. Masses that 
satisfy this condition are small in the sense that the resonance is shut off by the variance 
rather than the mass term in the lagrangian, and the maximum variances reached are hence 
the same as in the massless case. If, on the other hand, the resonance is shut off by the x 
mass before the scattering regime sets in, the maximum variances reached are smaller by a 
factor exp[— 2/iCc;°(rs — t^)], where Tm, the conformal time at which the mass term kills the 

resonance, is defined by = gy^((X'°)^(ao/a(rm))^. The borderline case is the maximum x 
mass for which the scattering regime is still reached: 

~ gQ a;^— ^ — - , (massless inflaton) . (32) 

This can be rewritten as an estimate for the minimum initial value of q required to reach 
the scattering regime for a given value of m^. Using 

In ^ , (33) 



^0 



and 

^=l + /f„a„r. ^„«_L_,;^, (34) 

where Hq is the Hubble parameter at the beginning of the oscillatory stage, one obtains 

min /^!%iLV ( ^ "^X "^0 1 1 ^scatt V /^c^ 

where Hq ~ 1/2, and n^catt ~ 9^^ for small A,^ and ~ 10/A^ for large A^ For go > Qo"^^ 
the variances peak at their scattering regime values. For a typical choice of parameters 
%/Mp = 1.8, = 0.1, and taking A^^ = 0, Eq. (^) yields q^"' ~ 2 x 10*^ for = cj°, and 
g™™ ~ 8 X 10^ for = lOu;^ |3^. This is roughly in agreement with Ref. [0, where the 
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authors estimate that no significant resonant production occurs for go < 10^ if < 00^. Of 
course the choice of /i introduces significant uncertainty into the estimates: using /i = 0.2, for 
example, lowers the values of g™™ by a factor of 5 — 10. To test the validity of Eq. ( pS] ) we ran 
our code with = l.OTtu^, = 0, and go = 3.65 x 10^ and found that the scattering regime 
is marginally reached when the resonance shuts off at t ~ 1.5 x lO^Mp ^. The peak and valley 
values of the x variance are ((5x)^)pcak ~ 4 x 10"-^°M| and ((5x)^)vaiiey ~ 6 x 10~-^^M|. This 
is in good agreement with the prediction for the scattering regime peak and valley values 
of Eqs. dll) and (P): ((5x)')peak ~ 4 x IO-^OM^ and ((5x)')vaiiey ~ 1-3 x lO-^M^, where 
we used for $s = ^o{ao/a{Ts)) ~ O.OSMp. Note that even though these variances are rather 
small the scattering regime is reached since large g corresponds to large g. For go = 10^ we 
find that the resonant growth shuts off well before the scattering regime is reached. The x 
variance grows to about 1% of its scattering regime values (|13[) and (|27|). 



2. Massive inflaton, g > 

Before discussing the effect of the x fi^ld mass on the resonance shut off we will discuss 
briefly the simpler situation with = but massive inflaton. Since for a massive inflaton 
$ oc a~^/^, we find that g = go(ao/Q-(^))^- There are then two possible mechanisms for 
terminating the resonance: either the resonance is shut off as usual by the variances, or g 
becomes <C 1 due to the expansion before this happens. In the former case the maximum 
variances can be estimated as usual from Eq. (|^). Hence Eqs. (|12D and (|13D for "small" A^^ 
(rg^/^ < (3c^)~^) and Eq. ([15| ) for "large" A^ apply also in the massive inflaton case if the 
variances shut off the resonance, i.e. if the scattering regime is reached ||3^. The difference 
is that $s is not obtained from Eq. (^2]) . Instead one must use 



a(ts) 



3/2 

(massive inflaton) (36) 



where 



a{t) f, ^rr V TT 1 "^0 1 , ^: 



l + --H,t]\ Ho^;-i^m^, t,^^ln%l, (37) 



Oq V 2 / 2 Mp 2m^fi n, 



and ts is the physical time at the beginning of the scattering regime. The estimates for the 
n^catt s-iid no are given below Eq. (pS]). 

On the other hand, if g <^ 1 before the variances reach the values obtained from the 
above estimates then the resonance shuts off simply because there is no strong instability 
for small g (c/. figure 1(a)). The regime with g ^ 1 is often called "narrow resonance" in 
the literature because the width of the instability bands goes to zero as g becomes small. 
Using the fact that g oc and the estimate for t^ above it is easy to find the minimum 
initial value of g such that g(ts) ~ 1- One obtains 

• K-O). (38) 

For typical values $0 = 1.8Mp and /i ~ 0.1 this evaluates to g™™ ~ lO'' in the small A^ 
case and g™™ ~ 461n^(|^) in the large A^^ case. If go ^ g™™ then only a small fraction 
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of the inflaton energy decays before q{t) <^ 1 and the fast decay shuts off because ^ 1. 
However, if go ~ •y™" then the inflaton can decay almost completely before entering either 
the scattering regime or the narrow resonance regime. This is the situation for which the 
variances of the fields reach their maximum possible values. 

When the x field is massive, one requires larger values of go for the variances to reach their 



scattering regime values. The minimum value of go for which the estimates of section II B 



together with Eqs. ( |5BD and (P7[), apply can be found as in the previous section by requiring 



I < q'^\Qml = q^\aQ/a{U)f'^ml. The result is 



If go is larger than this, the resonance is terminated by the variance. Taking the same 
parameter values as in the massless inflaton case of section |III B 1| we obtain g™™ ~ 10^ 
for = and g™™ ~ 3 x 10^ for = lOm^. Hence we conclude that for a massive 
inflaton the values of g™" are typically one to two orders of magnitude smaller than in 

1/2 

the massless inflaton case. A comparison of Eqs. and ( |5UD leads to ^massive/o'massiess 



{9c'^/2){'m^/uj^Y, where we assumed equal occupation numbers, /x's, x masses, $o's, and 
a;° = m^. We should point out that for a massive inflaton the slow roll ends somewhat 
earlier and that hence the appropriate values of $o are perhaps 20% larger than in the 
massless case. To summarize, the main reason why the massive inflaton is more efficient in 
producing massive particles is that the typical resonant momenta redshift slower. This gives 
the resonance more time to create particles before it is shut off. 

It is interesting to compare our analytical formula for the maximum variance reached 



with the estimate given in ||T3[. As explained above, the expressions given in section [II B 
apply as long as go is large enough so that the scattering regime is reached. For small 
the appropriate formula is then Eq. ([T3|), which is particularly easy to apply for a massive 
inflaton since the combination $(t)^/g(t) = Arn^/g does not redshift due to the expansion. 
Hence ((5x)^)peak ~ ^o/^QQo- For go = 10®, which as shown above is sufficient to reach 
the scattering regime even when m^/m^ = 10, we thus obtain {{SxY)peak ~ 10~^M|. This 
is in good agreement with [0, where the authors find ((5x)^)pcak ~ 10^^M| [|^. Note, 



however, that the parametric dependence of the estimate for {{6x)^)pea,k given in |]T3 
to be different from our Eq. 



seems 



3. Negative Coupling 

The case of inflaton decay into massive x particles with a negative cross-coupling g is 
signiflcantly more complex than the positive g case. We hence begin by discussing the main 
characteristics of the process. 

The flrst notable feature is that the resonance is shut off by the x mass only when 

ml > \g\^^ ^ A\q\cul , (40) 

which should be compared to the positive g criterion ^ q^^'^cu'^. Eq. (^) is simply 
Eq. (|10D for A{k = 0) = -2|g| + m^/cj? and 5A/A < 1, g > 1. Hence in order to get 
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significant resonant production, one requires <^ Ifl'I^O) which is a much weaker condition 
than for positive g. Since $ scales with the expansion of the Universe (as oc for a 
massless inflaton and oc a~^/^ for a massive inflaton), there is then a time interval during 
which resonant production is possible. The growth mechanism is quite different from the 
positive g case, and we will describe its main features below. However, a complete analytical 
understanding is extremely difficult to achieve. The reason is the extremely "chaotic" nature 
of the process, which we have already pointed out when describing the figures in section P"ll A . 

It turns out that the positive g resonant scale kf^/^ = |g|4u;^ plays an important role in 
understanding the dynamics. While <C \q\^^'^uj'^, the dynamics is equivalent to a positive 
g case with gcfr = 2|g| and /icg ~ 0.1, essentially because Xo evolves according to Eq. (||) 
(with ^ 0). Once > kf^/^, which always occurs sooner or later since kf^^ redshifts 
with the expansion of the Universe, the 'positive g resonance' shuts off completely. This 
does not mean that production stops altogether. There is another effect which we refer to 
as the negative coupling instability. 

We now give a heuristic explanation of this new effect. In some sense we already discussed 
it in section [II B| , where we argued that once the variances get so large that Srriy^ ~ \q\^^^ujfj,, 
Xo falls out of phase with xq ( cf. also [0). The difference is that for massive x ths "falling 
out of phase" can occur without the need for large variances. The reason is the following. 
After 00 passes through zero and starts growing, it tries to pull Xo along, attempting to 
keep it equal to Xo prescribed in Eq. (H). However, if is large enough so that Xo has 
time to oscillate while 0o is small, then a phase mismatch can occur and prevent xo from 
growing and following xo- This will be illustrated with a simple model below. When such 
a mismatch occurs the infrared x niodes grow very fast since for these modes the frequency 
squared is negative, and their potential corresponds to an inverted harmonic oscillator. The 
growth exponent /ieg = Xk^dxk/dt may be as large as ~ \q\^/'^iu^, and the range of momenta 
which grow can be as broad as Afcpj^yg ~ Is'l'^'^- A dramatic effect indeed, compared to the 
usual positive coupling parametric resonance. This is what we mean by negative coupling 
instability. The explosive growth ends once Xo grows sufficiently large to catch up with xq, 
at which point the infrared modes resume their normal oscillatory behavior. An example of 
this effect is presented in figure 11(a). For = 5.4 x 10~^^M| the variance grows with an 
exponent of order fx^s ~ 14 ^ 2\q\^^'^mfj,. 

A quantity which is clearly of interest is the probability that for a given set of parameters 
a negative coupling instability, and the huge particle production associated with it, occurs. 
As explained above, the instability can occur if there is a phase mismatch between xo and 
00 at the moment when xo should stop oscillating about the origin according to Eq. (||). 
Specifically, the instability can occur at the instant when the curvature of the x potential 
V{x) at the origin changes sign, which happens when 0o becomes larger than m'^/\g\. This 
observation can be used to estimate the a priori probability for an instability to occur, based 
on the assumption that the phase of xo is random at the instant the curvature changes sign. 
We stress that this assumption is reasonable only when m^gg ^ \q\^^'^uj'^^Q. If this condition 
is not satisfied Xo does not have time to oscillate while 0o is small, and hence no phase 
mismatch can occur. In this case Xo simply follows Xq, just as in the first part of the 
evolution shown in figure 2. 

One might guess that the a priori probability for the mode amplitudes to grow by a 
factor G ought to be of order G~^. This estimate is supported by the following toy model. 
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Let us take the time at which the curvature of V{x) at the origin changes sign to be t = 0. 
Assume that, while 0o — 0, X = xq satisfies 



X = Xq cos rriyt H sin m^t , for t < . 



(41) 



When 00 becomes larger than m^/\g\ at t = 0, Xo finds itself, at least momentarily, in the 
potential of an upside down oscillator. Hence we have 



^0 \ -ut , 1 



X, 



\ „vt 



for t > 



(42) 



where v = fitOfj, ~ We are interested in the possibility of a phase mismatch, by 

which we mean that X matches mostly onto the decaying solution and hence the amplitude 
of the growing mode at t = is small, i.e. 




eAr 



where 




e< 1 



(43) 



To recover its original amplitude Aq, and catch up with xq, X has to grow by a factor 
exp z/t ~ e^^ ^ 1. During this time period the infrared modes are unstable, and the x 
variance may grow by a factor as large as ~ ~ G^. The question we want to answer is, 
how likely is such an event? As mentioned above, the estimate we are about to present is 
based on the assumption that xo has a uniform random phase at t = 0, which should not 



be a bad approximation for > 



The appropriate probability distribution is then 



the uniform distribution on the circle of phases: 
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Xo = Aq cos (y9 , — ^ = Aq sin . 



(44) 



After some simple algebra we obtain the desired probability 



1 



2Ao 



Xo + ^ 

V 



< e 



4 



TT 



1 + {m^jvy 



e< 1 



(45) 



The probability that the variance grows by a factor ((5x)^)/((5x)^)o ~ is then ~ 
(2/7r)G~^ ~ C"^, as anticipated above. Here we used G ~ for the amplitude growth 



factor and the fact that production becomes impossible when vn? ^ Ifi'l'^*^ (c/. Eq. 



so Eq. (^51) apphes only when v 



\q\ ' ^<f> 



> 



m 



= 5.4xlO-"M| 



In figure 11(a) we present a run with 



75m0 in which the x variance grows by a factor 3 x 10^ in the interval 
t e [2.21,3.16] X lO^Mp^ According to our estimate the probability for such a large growth 
to occur is about 0.8%. Figure 11(b) shows xo vs. time for three slightly different x masses 
to illustrate the dramatic effect of a phase mismatch at the instant when the curvature of 
the X potential at the origin changes sign. The ml = 5.4 x 10~^^Mp run matches almost 
entirely on the decaying mode, so in that case Xo spends a long time near the origin and 
the variance grows by a huge factor. As can be seen in figures 11(a) and (b), a tiny (1%) 
change in the x mass causes a different matching and changes the resulting variances by 
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several orders of magnitude This explains the spiky nature of the curves in figures 7 
and 10. The high sensitivity of the growth to small changes in the parameters can be also 
seen in figure 12. Note that because of the chaotic nature of the production mechanism one 
should not take the fact that we have connected our data points in these figures with line 
segments too seriously - for all but the most closely spaced points there could be numerous 
peaks and valleys in between. 

The probability Ve in Eq. (^) is essentially the probability per inflaton half-period that 
the variance grows by a factor e~^, provided that q^^^io^es ~ "^xeff ~ Q^^'^^c^es- Note 
that, since the x variance contributes to rn^^f^, the last inequality also determines the 
maximum growth that can be achieved. This will be discussed in detail below. Here we 
want to emphasize that the probability is per inflaton half-period because a phase mismatch 
can occur twice per period when the curvature of V{x) at the origin changes sign from 
positive to negative. This means that to obtain the total probability for a negative coupling 
instability to occur for a given set of parameters one must multiply by twice the number 
of inflaton oscillations that take place while g^/^cu^efr ~ f^x^s ~ ^^^'^^(pes is satisfied 



Since we are mostly interested in very heavy particles that could not be produced with 
g > 0, we will assume that the first part of the inequality is satisfied from the start. It is 
then easy to calculate how large l^^l needs to be in order to have n chances for instability 
before g^/^ < my^es/^cpes and all production shuts off. The minimum \g\ for n chances is 

where /? ~ 37r/4 for a massive inflaton and P ~ 7r/2-\/2c for a massless inflaton. Since $o is 
also slightly larger in the massive case we see that one requires somewhat smaller values of 
|5f| in the massless case. 

We will now briefly discuss the ratio of "peak" to "valley" variances in the scattering 
regime for massive particles. This was discussed for massless particles at the end of sec- 



tion |1I B| , where we found distinctly different behaviors for positive and negative g. If the x 
field is massive the situation is more complex. For g > the ratio ((5x)^)pcak/(('^x)^)vaiiey 
is still given by ujmux/ ^^mm ~ Q^^^- Of course this is only valid if the scattering regime is 
reached. For g < two kinds of behaviors are possible. If the resonance is shut off by the 
X mass, i.e. the scattering regime is not reached, then the curvature at the origin of the 
X potential remains positive after resonance shut off and the oscillating 0o(^) simply mod- 
ulates the frequency of the x oscillations. In this case the frequency of the x field changes 
continuously, so the "random phase" approximation discussed at the end of section |11 Ij| ap- 



plies, and the ratio is again ((5x)^)peak/ ((^x)^)vaiioy ~ ^max/^min- The difference is that now 

^max = ^l + 3A;^((5x)^)valley ~ ^Tl^' ^min = "^x ~ 1^1"^^ + (('5x)^)peak ^ - |^|<l>2. 

If, on the other hand, the scattering regime is reached when ^ then the infrared 

X modes can be dragged along by 0o just as in the massless negative g case, and hence 

peak / ^ max valley 



((f^X)^) valley V t^min / ^ 



,2 



(47) 



which reduces to ~ {m^/uj^Y when > 3A;^(((5x)^)vaiiey, and to ~ \g\^^'^^/ujfi, when < 
3Ax((5x)^)vaiiey Here we used 3A^((5x)^)pcak ~ Ifl'l^^- Note that eventually the expansion of 
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the universe always reduces $ sufficiently so that the condition < |(7|$^ becomes violated 
and this case turns into the one discussed above Eq. (|47|) . 

We are now ready to discuss the specifics of the massless and massive infiaton cases. 



4. Massless Infiaton, g < 

In figure 7(a) we show the peak variances reached as a function of for various values 

of I go I- Figure 7(b) is a blow-up of the |go| = 350 curve. For these runs uj'^ = cV^</,$o ~ 
1.5 X 10~^Mp, so the horizontal range in figure 7(a) corresponds to m^/uj^ 0.63 — 66. The 
range in figure 7(b) corresponds to m^juP^ ~ 1.84 — 8. As discussed in section [111 B 1| , for 
positive g the scattering regime would never be reached for these mass ratios with go = 350. 
Even more dramatic is the fact that for m^/Mp ^ 2 x 10~^^, there would be no production 
at all for the chosen parameters if g were positive This illustrates in concrete terms the 
power of negative coupling to produce massive particles. 

The spiky nature of the curves was already explained in the previous section. In fig- 
ures 8(a) and (b) we show the variances and occupation numbers [B^] of a run with a x fi^ld 
mass for which there happen to be no phase mismatches and hence no negative coupling in- 
stability. In figures 9(a) and (b) we show the same quantities for a case where the instability 
does occur. Note that is actually larger for the run in figure 9. 

The minimum value of |go| for which there is any chance of production can be immediately 
obtained from Eq. (M) with n = \: 




Cn=ki 1 + • (48) 



With $0 = 1.8Mp this evaluates to Ig^'""] = 2.8 for m^/u;^ = 1 and Ig^'""] = 280 for 
m^/uj^ = 10. These values are orders of magnitude less than the corresponding positive g 



values derived in section [II B 1 



The maximum variance which can be reached is determined by Eq. (^) (with 

ml gfj) to be 

((^X)^)max = '"'^3""' . (49) 

The appropriate value of $ depends on exactly when the growth occurs. If a phase mismatch 
and the associated negative coupling instability occurs very early on (right after the first 
zero crossing, say) $ can be as large as 0.5Mp. 



5. Massive infiaton, g < 

When discussing the massive infiaton case with negative coupling it is important to keep 
in mind that the stability bound r > 1 always requires the presence of a X^(f)'^ (or higher 
order) term in the effective potential. For the runs in this section we have chosen the infiaton 
mass such that ^ 3A0$q at the end of slow roll. This clearly corresponds to mild fine 
tunning as there is no a priori reason that these quantities should be of the same order. 
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If one chooses <^ m'^/Mp it becomes more difficult to satisfy r > 1 and the available 
parameter space is somewhat limited. For these reasons we believe that with negative g 
the more natural situation is the one with -C 3A(^$05 which case the anlaysis of the 
previous section applies. Nevertheless it is interesting to study the situation where the 
inflaton dynamics is dominated by the mass term during the oscillatory stage, and we do so 
below. 

Figure 10 shows the maximum peak and valley variances reached during preheating 
vs. for I go I ~ 350 (the exact parameters are given in the caption). The horizontal 
range corresponds to 2.4 < m^/m^ < 12. According to the analysis of section p^llB 2| the 
scattering regime is never reached for the same parameter values with g > 0. Moreover, for 



m^/Mp ^ 10^^^ there is no production at all for the chosen value of |ff| if > |^ . 

As in the massless inflaton case the minimum value of |go| for which there is any chance 
of heavy particle production is obtained from Eq. ( ^61) with n = 1: 




IC1 = hE^ 1 + ^^ • (50) 



For m-^/m^ = 1 this evaluates to Ig™"] = 6.9, and for m-^/m^ = 10 one obtains |g™™| = 690 
(here we have taken $o = 1.8Mp). Comparison with the results of section [III B 2| shows that 
once again the required values are orders of magnitude less for negative coupling than for 



positive coupling. On the other hand comparison with section [1IIB4| shows that for negative 
coupling the massless inflaton is actually somewhat more effective in producing massive 
particles than the massive inflaton. This is opposite to the positive coupling situation, 
although it should be noted that for the present case the difference arises only from a 
numerical factor; the parametric dependence is the same. Finally we point out that the 
maximum possible variance is given by Eq. (|4^) also for the massive inflaton since Eq. (^OD 
holds in both cases. 



IV. BARYOGENESIS 

It is well known that grand unified theories generically predict the existence of baryon 
number (B) violating interactions. Much effort has gone into converting this fact into a viable 
mechanism for baryogenesis In the standard scenario, B violating gauge or Higgs bosons 
fall out of thermal equilibrium as the temperature of the universe drops below their mass 
and subsequently decay, creating net B in the process. This B must then be protected from 
being wiped out by sphalerons, electroweak B violating processes. Fortunately sphaleron 
transitions conserve B-L, where L is lepton number, so that any GUT that creates a non 
vanishing B-L remains a viable candidate. Examples are 50(10) and Eq. 

According to the old theory of reheating , based on perturbative decay of the inflaton, 
GUT baryogenesis seems rather difficult to implement [^. The reason is simply that the 



heavy B violating gauge or Higgs bosons cannot be produced efficiently. The gauge bosons 
are expected to have masses of order the GUT scale (lO^^GeV), while the Higgs bosons may 
have masses as low as lO^^'GeV. Since the mass of the inflaton cannot be larger than 
about lO^^GeV in order not to be in conflict with the measured anisotropy in the microwave 



background direct perturbative production of GUT bosons is kinematically forbidden. 
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In addition, the weak coupling of the inflaton to other fields causes the reheat temperature, 
as obtained in perturbation theory, to be so low that thermal production is also ineffective 

Clearly GUT baryogenesis after infiation is resurrected by parametric resonance. Not 
only can bosons with masses larger than the infiaton be produced directly. They are auto- 
matically far from thermal equilibrium, thus satisfying one of the prerequisites for baryogen- 
esis. The first model following the baryon asymmetry from the end of preheating to the final 
decay of the infiaton is presented in |^8[. The authors assume that during the resonance 



stage a fraction of the initial infiaton energy is rapidly transferred to GUT bosons with mass 
lOm^. These bosons then decay in a CP violating manner, producing net baryon number 
in the process. 

While the results of the toy model in are very encouraging, it has a rather unpleasant 
feature. As shown in as well as sections [III B 1| and P^IIB2| , for (7 > it is quite difficult 
to produce particles as massive as lOm^ during preheating. Specifically, production of such 
particles requires a massive infiaton and go ~ 10^. Since ^ 10^'^GeV this corresponds 
to g ^ 10~^, which leads to a fine tunning problem: in order to reproduce the anisotropies 
in the microwave background the quartic self coupling of the infiaton is constrained to be 

^ 10~^^, which is unnatural because A<^ receives contributions of order g"^ at the one-loop 
level. As mentioned in |H8| this problem does not arise if the value of A,^ is protected by 



a supersymmetric non-renormalization theorem. However, so far no such model has been 
constructed, and in any case it is interesting to explore other solutions to the naturalness 
problem. 

In light of these considerations our investigation of massive particle production in sec- 
tion Q] is particularly important for baryogenesis. We have shown that with g < large 
amounts of heavy particles are produced even for moderate values of go- In particular the 
naturalness condition g'^ ^ A,^ is easily satisfied even for ~ lO^^GeV or more. Addition- 
ally, the production of the heavy particles is extremely rapid, so that it is not affected by 
competing decay channels. In fact, as explained in section |III B| , for g < the decay into 
massive particles is much faster than the decay into light particles provided the mass is such 
that xo and Xo quickly fall out of phase. Finally, in contrast to the positive g case, with 
g < heavy particles are produced by both massive and massless inflaton flelds. For these 
reasons we feel that negative coupling baryogenesis after preheating is an exciting possibility. 

In order to illustrate how this baryogenesis mechanism works we adopt the toy model 
presented in [Q, with some minor modiflcations. Before doing so let us briefly mention 
two potential problems with the GUT baryogenesis scenario in general. First, there is the 
possibility that GUT symmetry is restored after inflation by a thermal or nonthermal |TB,n] 
phase transition. This could lead to the production of cosmologically dangerous topological 
defects. Fortunately it seems unlikely that either type of phase transition would be strong 
enough to restore the symmetry, at least for large portions of parameter space |I9| , ^ . The 
second potential problem with the scenario is that gravitinos are overproduced and spoil 
nucleosynthesis. This also does not look like a major obstacle because the late decay of 
the inflaton has the ability to sufficiently dilute the gravitino abundance produced early on 

mm- 

The toy baryogenesis model we will analyze is as follows. When the production of massive 
particles shuts off ($ = M^/ ^J\g\) [^, there is a fraction of energy 5o = p% in the GUT 
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boson field x- As indicated by our results in section |T|, we will take 5q <^ 1. Since the x 
bosons are only moderately relativistic at creation {k ~ M^) they quickly become "cold" as 
the Universe expands and can be approximated as having zero kinetic energy. We assume 
that they decay into light particles with a constant decay rate F^. The decay of each x~ X 
pair produces e baryons, where e should be thought of as an effective CP violating parameter 
4^ . In a reahstic model e is expected to be small. We also assume that the decay products 



thermalize instantly, which is a good approximation if they are strongly coupled. The decay 
of what remains of the infiaton field is modeled by a constant rate F^, and we will assume 
a hierarchy of time scales 

H,^ « F-i « F^i , (51) 
where Hq is the initial Hubble parameter. Finally, we make one important simplification 



over the model presented in |^8[: it is assumed that at all times the temperature of the 
decay products is much smaller than the mass M^. This allows us to compute the baryon 
to entropy ratio ub/s analytically. The precise criterion for the validity of the inequality 
^ Tmax will be given below. 

With these assumptions, the x ^'^^ baryon number densities per comoving volume, = 



n 



[a/aoY and A''^ = nsia/aoY, evolve according to 



4 = -F^iV, , Nb = -eiV, , (52) 

with initial conditions NB{t = 0) = and N-^{t = 0) = A^^ = n°. Note that we have 
neglected the production of x particles via inverse decays, which is justified because of 
our assumption that the temperature T -C at all times. This condition will be made 
quantitative below. The solutions for Nb and A'^ are 

= N^e-^^' , Nb = eiV° (l - e"^^*) (53) 

The first step in computing the baryon to entropy ratio is to solve for the energy density of the 
baryon fluid, Prad, which is related to the entropy density s by s = {2tt'^ / A5){30 /n'^y^'^ gl^^ p'^[^, 
since the baryon fluid is relativistic and thermal. Here g^ is the number of relativistic degrees 
of freedom, and we are assuming that the entropy is dominated by the relativistic particles. 
We will consider separately the massless inflaton, whose energy density scales as oc a~^, 
and the massive infiaton, with oc a~^. 

For the massless infiaton we will assume that the universe is radiation dominated at all 
times, which requires 3> p° (i^o/r^)^^^. Then 



P.ad f ) = -^M,N, - I L 




a 



1 4> 



ao 



(54) 



= {2Hot + l)K ^o = ^p5- (55) 

Note that in order to compute the late time contribution of the last term on the right 
hand side of (p4|) to prad, we need not know the details of how the infiaton decays, i.e. 
the exact time dependence of p^^a/a^Y. The relevant information is the time scale F7^ at 
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which the decay takes place. We will first study the solution to Eq. ( p4D at early times 
{t -C -C r|^^) to obtain an expression for the maximum temperature. This will allow us 
to derive the inequality that needs to be satisfied in order that T^ax M^. Subsequently 
we will obtain the solution at intermediate (r~^ <t< F^^) and late times {t > T^^), which 
will be used to compute ub/s. 

Early times [t -C T^^) are relevant since, as we will see in a moment, the temperature 
reaches its maximum value at t ~ H^^ <^ Fr;^. Using the fact that in our simple model 



(ao/o)'^P0 exp(— r^t) and early on 



-r^A^o, we find 



Prad 



{a/aof 



Ho 3 (a/ao)- 



1 {a/aof 



{a/aoY 



(56) 



Note that for a/ao ^ 1 the two terms on the right hand side scale as a ^ and a ^, respec- 



tively. Hence pmax is reached when a/ao is of order one, i.e. t 



Ho \ 



Pmax — 30^* 



4-^/Vy§^ , at a = 4^ , when p^V^ > 




2t 



X 

when Px^x < 



and we obtain 



(57) 



so that 



< Mx if 



(58) 



In this paper we assume that Eq. (]58|) is satisfied. In addition to Eq. (pSj) there is an 
upper bound on the baryon to entropy ratio ub/s which can be used as a consistency check. 
Since we assume instant thermalization of the baryon fluid and Tmax < M^, the number 
density of relativistic particles g*nra,d is always larger than the number density of decayed 
X-X pairs, n^/e. Since the number density of relativistic particles is related to their entropy 
via rZi-ad = [C(2)/4C(4)]s ~ 0.28s (where C,{z) is the Riemann zeta function) we hence have 



s 



< e 



C(3) 
4C(4) • 



(59) 



We will now derive an expression for the baryon to entropy ratio in our model. At interme- 
diate (r^^ > t > r~^) and late times {t > F^^) the scale factor is to a good approximation 
a ~ {2Hotf''^. Hence from Eq. (||) 



a 

Prad \ 




exp 



F. 



2Ha 



a y 



\0'Qy 



(60) 



where we used /o°°(itt^/^exp(— t) = 7r^/^/2. The second term on the right hand side of 
Eq. (^) first grows oc and then, when the inflaton decays at t ~ F^^, it attains a 
constant value p°. Since we are assuming that p° ^ p° (i/o/rx)^'^^ ^"^^ "final" 
baryon to entropy ratio 



Ub 




(massless inflaton) , 



(61) 
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where we used A^^ = for t ^ F"^. For simphcity in this, and in the results below, we 
ignore multiplicative numerical constants that change the final result by less than about 20%. 
The remarkable feature of Eq. ( ]6T| ) is that the final baryon to entropy ratio does not depend 
on the decay constants F^ and F^. Why this is so can be understood as follows. As long as the 
temperature never rises close to M^, so that the massive x particles are never re-populated 
by inverse decays, they simply decay out of equilibrium and create Nb — €N^ baryons, 
which is independent of the decay rate F^. The entropy produced per comoving volume, if 
dominated by the infiaton decay, is s(a/ao)^ — [^^(a/ao)^]^'^^ — (p^)^^'^, again independent 
of the decay rate F,^. Hence the final baryon to entropy ratio is also independent of the decay 
rates. In particular, even if the decay is not as simple as the one-particle out-of-equilibrium 
decay, but a more complex process like parametric resonance, we expect Eq. (|6TD to be valid 
as long as the infiaton decay time is much greater than the massive particles decay time, 
F^ <^ F^, and the entropy production is dominated by the late time infiaton decay. This 
concludes our discussion of the massless infiaton, and we now turn to the massive case. 

When the infiaton is massive, its number density per comoving volume decays as A^^ = 
ntpi^a/aoY = A^° exp(— F^t). The evolution equation for the radiation density is then 




" a ' 




a / 


.ao. 




ao ^ 



(62) 



Note that because of our assumption that p° ^ p° the universe will remain matter domi- 
nated until the infiaton decays at t ~ F^^"*^. At early times (t <^ <^ F_^^) ^ — F^A^°, 
and hence 

Prad 



^0 Ix , „o £0 



5 (a/ao)^ ■ ^ ' 



This peaks at a/ao = (8/3)^/^ ~ 1.48 (t ~ 0.52ifQ ^), leading to a maximum temperature 



T ( T 

^ max / rx X 



Thus Tinax ^ requires 



plT^ + plT^^g^'^Ho, (66) 



which is the same constraint as in the massless infiaton case, Eq. (|58D. 

At intermediate times (F^^ > t > F~^) the scale factor is approximately a/ao 
{2>Hot/2f''^ and hence 




where we used dtt'^^^ exp{—t) = (2/3)F(2/3) ^ 0.90. This leads to a baryon to entropy 
ratio 
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While the second term in square brackets is subdominant, n^/s is to a good approximation 
constant and determined by the entropy release from the massive particle decay. Eventually 
the second term takes over and ub/s oc a~^^^^ until the inflaton decays at t ~ F^^ 
At late times {t > F^^) the radiation energy density is 

p-f-) =(^) r(^)(^) pI (69) 



leading to a "final" baryon to entropy ratio 

Ub 





(massive inflaton). (70) 

This means that for a massive inflaton the final value of the baryon to entropy ratio is 
suppressed by a factor (F^/ifo)^^^ in comparison to the massless case, Eq. (|61|). Note 
also that ub/s is dependent on the decay rate F<^, unlike in the massless case, but still 
independent of F^. Intuitively, this dependence can be understood as follows. For a massive 
inflaton the energy density scales as a~^, which is slower than for the massless baryonic fluid 
(oc a~^). Hence, later inflaton decay leads to more abundant entropy production relative to 
the entropy of the baryonic fluid. 

Before giving an estimate of n^/s based on the equations derived above we would like 
to emphasize that the assumptions Tmax ^ and p° 3> p° are in no way necessary 
for successful GUT baryogenesis after preheating. These conditions were imposed only to 
enable us to obtain simple analytical expressions for the baryon to entropy ratio. For a more 
sophisticated numerical treatment in which both these assumptions are relaxed see . 



In order to obtain estimates for the baryon to entropy ratios that could be produced at 
preheating in a theory with negative coupling we need to know p^/ p^ at the time when the 
massive particle production shuts off. To this end we ran our code for a range of "natural" 
g values with mx = lO^^GeV. The results are shown in figures 12(a) for a massless inflaton 
and 12(b) for a massive inflaton. We used the COBE values = 3 x 10^^'^ in the former 
case and = 3 x 10^^'^, = lO^^GeV in the latter case. We point out that for natural 
values of g there is no production with these parameters if g > 0. The coupling A^. was 
chosen so that r = 10 in all of the runs. 

The ratio of energy densities p°/p^ can be easily obtained from figure 12. To see this 
recall that for the massive nonrelativistic x particles the variance is simply related to the 
energy density by p^ = M^((5x)^). Now if there is any growth at all we find that ((^x)^) does 
not decrease until the time at which the massive particle production is shut off. Hence the 
values in figure 12 are proportional to p° , the energy density at the instant when production 
ceases. The corresponding value can be obtained by estimating $ at this time from 
Eq. (|40|). Note that this last step is a good approximation only for very heavy particles such 
as the ones we are considering here. For lighter x particles the 3A^((5x)^) contribution to 

eff dominates over when the production shuts off. 
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Here we will concentrate on the massless inflaton case and obtain an estimate for ub/ s 
based on Eq. (|6l|) Using = M^{{6x)'^)o and = A^($0)74 ^ X^M^/A\g\^ we obtain 

From figure 12(a) we see that the largest valley values are about ((5x)^)o/^p ~ 5 x 10""'^°. 
It is easy to check that for these values the constraint Eq. (|58| ) is indeed satisfied provided 
r^, <C Hq. Using this value for (((5x)^)o in Eq. ([71]) we find 

~ 10-3 , (72) 

e s 

where we have set g^, = 100. This leaves plenty of parameter space for baryogenesis consistent 
with the observed value (tt-b/s) observed ~ 4 — 7 x 10"^^. 



V. SUMMARY AND CONCLUDING REMARKS 

In this work we investigate the preheating dynamics in chaotic inflationary models with 
a negative cross-coupling between the inflaton and a second scalar field. Such couplings 
are completely natural in theories with more than one scalar multiplet. The details of the 
inflaton decay are extremely complicated, and an accurate treatment requires that one solve 
the full nonlinear problem. Fortunately the dynamics is dominated by states with large 
occupation numbers which admit a semiclassical description, so the problem can be treated 
by solving the classical equations of motion numerically on the lattice. We use a general 
renormalizable two field scalar theory that respects a Z2 x Z2 discrete symmetry, i.e. the 
potential contains quadratic and quartic terms as in Eq. (|l]). Perhaps our most important 
finding is that the negative coupling can catalyze the inflaton decay into a heavy scalar field, 
opening a new window for baryogenesis at the grand unified scale with a natural choice of 
parameters. 

In order to make our presentation self-contained, we review some of the important prop- 
erties of inflaton decay for massless fields with positive coupling in section |I|. In the old 
theory of reheating with perturbative inflaton decay the relevant quantity was assumed 



to be the reheat temperature. Based on the reheat temperature one can infer whether the 
grand unified symmetry is restored, whether there is a monopole problem, whether GUT 
baryogenesis is feasible, etc. In the new theory of reheating the inflaton decay products are 
typically low energy, far from equilibrium excitations, so that one cannot meaningfully as- 
sign a temperature to this infrared condensate. Instead, the relevant measure of the inflaton 
decay efficiency is the maximum variance the decay products reach since, with the variance 
known, all of the problems mentioned above can be addressed. Hence we focus on obtain- 
ing estimates of the maximum variances reached during preheating. The numerical results 
for cross-coupling g > are presented in figures 3 and 5, and estimates for the maximum 
variances are given in Eqs. ([T^ ) - ([T6|) . (Recall that is the inflaton and x is a scalar field 
coupled to it). The estimates are based on the observation that for much of the relevant pa- 
rameter space the backreaction of the produced particles terminates the exponential growth 
when the variances reach a certain magnitude, which can be understood intuitively in terms 
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of the stability chart of the Mathieu equation, shown in figures 1(a) and (b). Once the 
critical variance is reached, a slowly varying state, which we call the scattering regime, sets 
in. At this stage of the evolution the slow particle production cannot keep up with the ex- 
pansion of the universe and the variances start decreasing. As explained in section [II B| and 



r^, the fact that our nonlinear lattice calculations take scatterings into account is crucial 



for the determination of the correct maximum variances. An important consequence of the 
existence of the scattering regime for large values of g and A^,. is the fact that the inflaton 
decays more efficiently for smaller couplings. For example, a massless inflaton decays fastest 
and produces the largest variances if the cross-coupling satisfies < (7 < 10^ A<^ (provided 
the self coupling of the decay product x is not too large, A^^ < g^^"^ / ^yj"^)- The explanation 
of this rather counter intuitive result - that larger couplings lead to less efficient decay - 
is simply that larger couplings lead to more efficient scatterings and stronger backreaction, 
and consequently lower maximum variances. (For completeness we note that for a massive 
inflaton the most efficient decay into massless x particles occurs when go ~ ^ZcT™ given in 
Eq. (^), which translates into 10^^ ^ g ^ 10~^ in the case when A^ is small. Here we took 
for the inflaton mass ^ 2 x lO^^GeV.) 

We begin our investigation of the negative cross- coupling case by considering the instabil- 
ity chart of the Mathieu equation (see Eqs. (Q) and (^, and flgure 1), which is essentially an 
oscillator equation with a time dependent frequency. The Mathieu equation is a good start- 
ing point since it describes the simplest possible case, namely the linearized mode equations 
of the second fleld in a static space-time when the inflaton is massive. Its usefulness goes far 
beyond that simple case, however. The reason is that the main features of the stability chart 
are extremely robust and carry over to both other models (such as a massless inflaton) as 
well as the expanding universe. The Mathieu equation contains the two parameters A and q 
which parametrize the frequency squared of the oscillator. The third important parameter 
is the instability index fi. The meaning of fi is simply that in one oscillation of the inflaton, 
the amplitude of an unstable x mode increases by a factor e^'^'^. In flgures 1(a) and (b) 
we show the regions of stable and unstable solutions as a function of A and q, as well as 
some curves of constant /x. The only features of these charts which we need in order to 
obtain analytical estimates for the maximum variances are (i) that the maximal distance to 
the flrst instability band above the line A = 2\q\ is approximately ^g, and (ii) that /i is a 
rapidly decreasing function of A above the line A = 2|g|. As stated above, these features are 
extremely robust, and from them one can derive the resonance shut-off condition Eq. (p!0D 
on which all of our estimates, for both positive and negative coupling, are based. 

For a quartic interaction term (yf0^x^/2 with a positive cross-coupling g, the phase space 
is limited to A > 2\q\. As a consequence the instability exponent fi ^ 0.3 for any q, which 
can be seen nicely in flgures 1(a) and (b). Notice also that in the flrst unstable band above 
A = 2\q\ n may peak at a value between n ^ Q.l and /i ~ 0.3, depending on the exact 
value of |g|. On the other hand, for a negative cross-coupling, all of the parameter space 
above A = —2\q\ is accessible. This means that fx < (4/7r)|g|i/2 > 1 (c/. Eq. (g)), and 
that there is a wide range of unstable modes: Afcphysicai ~ 2|g|^/^c<j0. The possibility of an 
extremely rapid decay into a broad range of momenta motivated us to study the negative 
cross-coupling case in detail. 

The following are our main observations and flndings. As explained above, naively 
the maximum value of the instability coefficient is //max — (4/vr)|g|^/^. However, such an 
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instability occurs only if the expectation value of the second field is zero at the end of 
inflation, i.e. Xo = 0- realistic situation this is not the case since for g negative the 
true minimum of the x potential energy is displaced from the origin (note that in order 
to have a potential which is bounded from below the stability condition Eq. (P) must be 
satisfied). Our analytical analysis and numerical simulations show that quite generically 
the dynamics of the fields during inflation will drive xo towards its true minimum. This 
implies that the natural initial condition at the end of inflation is Xo = Xo? where xo is the 
(00 dependent) location of the minimum of the x potential, given by Eq. (0). Using this 
initial condition we studied the case of massless fields with negative coupling in section 
We found that during the early stages of the post inflationary evolution the unstable modes 
grow just as in the positive coupling case, with instability coefficient fi ~ 0.2. This can be 
observed in figures 2 and 4. The reason for this behavior is as follows. As long as the x 
zeromode follows its true minimum xo = Xo? no x niodes encounter genuine negative coupling 
instability (negative frequency squared), and hence they grow only via the standard positive 
g resonance mechanism, with q^g ~ 2|g|. This picture breaks down once the backreaction of 
the created x particles becomes large enough to disrupt the evolution of xo so that it stops 
following Xo- When this happens the mode amplitudes become unstable and grow very fast, 
quickly reaching values large enough to terminate the instability, and the slowly varying 
scattering regime sets in. During this stage we find that the dynamics of the negative cross- 
coupling model differs profoundly from that with a positive cross-coupling. Perhaps the most 
striking difference is that the x variance peaks at a significantly larger value - by a factor 
4|g|^/^ larger than in the positive g case - leading to a much more complete inflaton decay 
into X particles. This difference can be understood from the instability chart in figure 1. To 
shut off the resonance in the positive g case, one needs a backreaction of order SA ~ g^/^, 
while in the negative g case a much larger shift 6 A ~ 4|g| is required. Estimates for the 
maximum variances with a negative cross-coupling are given in Eqs. (^O]) and (0). Another 
notable difference between g > and g < is the way the variances fluctuate between 
"peak" and "valley" values as the inflaton oscillates in the scattering regime. We present a 
simple explanation of this phenomenon using the matching of two oscillator solutions with 
different frequencies (see Eqs. ( p^ - (pHD). The peak to valley ratio in the case of massive 
fields is discussed at the end of section [III B 3| . 

In section |T| we study in detail the decay of the inflaton into massive particles (m^ ^ 0). 
Sections [III B 1| and [III B 2[ discuss the positive coupling case for a massless and massive 
inflaton, respectively. The main conclusion is that the estimates for the maximum variances 
obtained in section [II B[ (Eqs. ([1^), (f3|), and ([15[)) hold also for a massive x field provided 
that go is larger than a critical value g™™ given by Eqs. ( P^D and (|5D|). (see also P5[].) 



For g < the situation is very different. The criterion for the complete particle produc- 
tion shut-down is > 1(71$^ (c/. Eq. ([40|) ), where $ = $(t) is the inflaton amplitude, which 
decreases due to the expansion of the Universe. As explained in section [IIIB 3[ , particles can 
be produced only at the moment when = \g\(f)Q{t), and while (po^t) is growing. The min- 
imum condition to encounter any massive particle production is hence \g\ > m^/$^, where 
$1 ~ lO^^GeV is the inflaton amplitude at the first extremum in the oscillatory regime. 
For = lO^^GeV, this gives \g\ ^ 10~^, which agrees well with our numerical results in 
figures 12(a) and (b). To encounter n chances to grow, the required value of l^^l increases as 
given in Eq. (HB|). Note that this value is many orders of magnitude smaller than the value 
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g ^ 10~^ required for the production of equally massive particles in the positive g case. The 
particle production mechanism is very different from the usual positive g resonance, and we 
call it negative coupling instability. The details are complicated, but the main mechanism 
can be readily understood and is explained in sections ^11 B 3| - |III B 5| . Characteristic fea- 
tures of this new process are that massive particles are produced in extremely rapid bursts 
over a broad range of momenta (much faster and broader than in the positive g case) . This 
is illustrated in figures 9 and 11(a). Another important feature are the sharp spikes in the 
production of massive x particles as a function of the parameters, seen in figures 7, 10, and 
12. This behavior is explained with the aid of a simple analytical model in section [111B3 



(see Eqs. (^) - (0)). The system may be considered "chaotic" in the sense that it is quite 
impossible to predict the exact position and amplitude of the spikes as a function of the 
parameters. 

The main result of section |T| can be re-stated as follows. For a realistic choice of 
couplings, which in chaotic inflationary models is constrained by the COBE satellite mea- 
surements to be ~ 3 X 10^^^ or ^ 2 x lO^^GeV, one can produce massive particles 
with ~ lO^^GeV (as required by GUT baryogenesis models), provided \g\ ^ 10"*^ and 

> 3 X 10~^ (for stability) (see figures 12(a) and (b)). This means that the parame- 
ter space available for production of massive particles with "natural" coupling constants 
is appreciable: 10~^ \g\ 10~^, 1 > > 3 x 10""^, leaving plenty of opportunity for 
baryogenesis model building. This is to be contrasted with the positive g case for which 
g ^ 10^^ is required. Such a large value leads to an unpleasant fine tuning problem since 
the small value of A<^ needs to be protected against radiative corrections. 

To obtain an estimate of the baryon asymmetry that could be produced during preheating 
we study a simple toy model in section In short, the model can be described as follows. 
Initially a certain amount of energy is transferred from the infiaton to a heavy GUT scalar 
field via the nonperturbative mechanisms described in this paper. The result is a cold (far 
from equilibrium) fluid of massive particles (with a mass M^), which we assume decays in a 
B and CP violating manner into light degrees of freedom that instantly thermalize. In order 
to treat the problem analytically we further assume that at all times the temperature T of 
the light relativistic particles is -C M^, which requires that Eq. ( |58D is satisfied. The main 
results of this investigation are the final baryon to entropy ratios given in Eqs. (0) and (|70D 
for a massless and massive infiaton, respectively. These equations are obtained assuming 
a "natural" hierarchy of time scales, i^o 3> S> F<^. As can be seen from Eq. (|6l|), 
under these conditions the baryon to entropy ratio for a massless infiaton has the pleasant 
property that it is independent of the decay rates F^ and F,^. We believe that this feature, 
which can be understood intuitively and is explained in section |V], lends more credibility 
to our result since the final estimate is somewhat independent of the details of the model. 
Eqs. (|6l|) and (^) say that the final baryon to entropy ratio is proportional to the initial 
energy density in x particles. Since the energy density of the (nonrelativistic) x particles is 
proportional to their variance (p^ ~ M^(((5x)^)), the spikes in figure 12 are directly imaged 
to spikes in baryon production. For a massless infiaton the final baryon to entropy ratio is 
expressed in terms of the maximum variance in Eq. (^TJ) . Our numerical results indicate 
that it is not hard to obtain values in the range ne/s ~ e x 10~^, where e is the effective 
CP violating parameter in our model. Clearly more realistic models for GUT baryogenesis 
at preheating will have to be constructed to verify the viabihty of this scenario, but the 
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preliminary results look promising. 

Another interesting question that deserves further study is whether it is possible to 
construct a model such that the negative coupling instability discussed in this paper does 
not require any "phase mismatch" to become operational. In other words, the aim is to 
obtain natural initial conditions at the end of inflation that immediately allow the infiaton 
to decay rapidly into a heavy scalar field. The preliminary investigation of one realization 
of such a scenario is underway. 

Close to completion of this work we received a preprint |[hep-ph / 9704452|] [p^ in which 
the preheating dynamics in an expanding universe is also studied in detail. In this very 
interesting paper the authors focus exclusively on the massive infiaton, positive g, small 
case, for which their conclusions partially overlap with ours. A detailed comparison will be 
deferred to future publications. 
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one obtains ^{rf) ~ 2c/i0$o/ ln(A;^^). For = lO'^^, % = l.SMp, and fi<p ~ 0.04 one 
finds {m')sc.tt ~ 10-6M|. 

[26] The fraction of the total energy in x fluctuations at the beginning of the scattering 
regime is of order 0.2 x for small A^, and of order 0.2 x g/X^ for large A^^. These 
estimates apply to both massive and massless particles. For a massless inflaton the last 
expression can be rewritten as ~ 0.1 x [rq)~^. In deriving these expressions, we used the 
fact that the energy density in fluctuations and the variance are at the beginning of the 
scattering regime dominated by the same momenta, and hence ~ (A;^es)^((^x)^)peak- 

[27] In our runs we observe that even after xo stops following xo its dynamics is still strongly 
influenced by 0o and it often oscillates with frequency nc<j<^/2, where typically n = 3,4 
or 5. For n odd, xo develops an expectation value, indicating symmetry breaking. Here 
we will not investigate further the physical implications of this observation. 

[28] For negative g, 3A^((5x)^) always dominates over \g\{{6(j)Y) for g ^ 1 due to the stability 
bound r > 1 (c/. Eq. (^). Here we are concerned with the early evolution before the 
X production shuts off. At late times can decay into its own fluctuations and hence 
((50)^) can grow signiflcantly larger than {{Sx)"^), but this is obviously irrelevant for 
obtaining an estimate of {{SxY)peak- 

[29] Equipartition of energy implies that pscf, ~ p^-^. The energies can be related to the 
variances as follows. In terms of the occupation numbers, the variances can be writ- 
ten as: (((5x)')scatt ~ i2nY'nUkAk)^ ~ (27r2)-in^(2|g|i/2)^2 ^ 

(kAk)^ ~ (27r^) '^nf.ool, where we have assumed that the variances are dom- 
inated by the momenta of order fc^ ~ Ak^ ~ fcres — {M'l\Y^^^<t> and k^ ~ Ak^p ~ to^, 
respectively. At early stages of the scattering regime the same momenta dominate the 
energy densities, so we can write ~ {2TT^)~^n^{k^Ak)^ ~ (27r^)~^n^(4|g|)tt;| ~ 

2\q\'/^ul{{6x)^) and ~ {2n^y'nt{eAk)^ ~ {27i^)-'ntuj^^ ~ This then 

leads immediately to ((50)^) ~ 2|g|^/^((5x)^), implying Eq. (^. 

[30] The maximum unstable momentum of the x fl^ld is approximately k^^^ ~ cu^g^^^ in the 
negative g case, while kf^^ ~ cu^. The ratio hence goes like g^/^, which requires excellent 
IR resolution in order to capture both accurately. This must be combined with sufficient 
ultraviolet (UV) range to accommodate the created particles and model their scattering 
towards higher momenta. Since the x particles scatter more efficiently for larger r the 
UV cutoff must be increased as r is increased. 

[31] Our numerical results indicate that the dependence of rigcatt on the parameters is actually 
more complicated, but the expressions given in the text are roughly correct. Since the 
dependence of Tg on nscatt is anyway only logarithmic and since as discussed in the 
text there is significant variation in the relevant value of p for small changes in the 
parameters the uncertainty in the value of rigcatt need not concern us here. 

[32] When gcj)^ is large, Xo is trapped in the deep minimum at Xo- However, as 0o approaches 
zero it looses its grip on xo and Xo wants to oscillate about the origin with frequency 
uj^ ~ {3X^{{Sxy)y^'^- In the beginning cu^ is very small and Xo does not have time to 
complete a significant portion of an oscillation in the interval when 0o is small (more 
precisely, when gcpl < 3X^{{Sxy))- It is easy to show that when 3X^{{5xy) ~ 2(7-^/^0;^$ 
the X zeromode has enough time to complete half an oscillation. It is then completely 
out of phase with xq when 00 becomes large again, and consequently stops following its 
minimum. Note that even before this complete phase mismatch occurs there are small 
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mismatches each time 0o grows after crossing zero, once (((^x)^) becomes appreciable. 
This is the reason why the x variance can grow larger than its positive g value given 
in Eq. ( ]T6| ) before Xq stops following Xo- We point out that the mechanism of Xo falling 
out of phase with xo will play an important role when we discuss the massive x case in 
section |ITT[ As we will see, in that case the effects can be even more dramatic. 

[33] As we will see in section |ITT| , this is no longer correct when the x fi^ld is massive. 

[34] This "random phase" approximation is clearly only sensible if there is enough time 
between matchings for the solution to "oscillate" . This requires that the n^^ time interval 
At„ between matchings satisfies At„ ^ l/t^n = e^'^/cUmax- The minimum total time 
required is then T = Y.n=o ^'^n = t^max[(^max/t^mm - l)/(e^ - !)]• In the text we are 
interested in the limit when e is small and the ratio of frequencies is large, so T ~ 
(ecUmin)"^- Now from Eq. ( pSj) we found 0;^;^ ~ g^/^cu^, and the total time for the 
frequency to change between cUmax and cUmin is a quarter period of 0o; T = n/iT.uj^). 
This yields e ~ Q'"^''^, and so the approximation in the text is reasonable for large q. 

[35] Note that we are applying the COBE constraints to the couplings of the field even- 
though the slow roll during inflation takes place in the direction (0o, Xo)- This is justified 
because the field dominates the dynamics during inflation, in the sense that the en- 
ergy density and the slope of the potential at the point (0O) Xo) are given by X^cj)^/^ and 
A(^0Q, respectively, up to corrections of order 1/r. Since r is larger than one (c/. Eq. (||)), 
and there is no reason for r to be close to this bound, these corrections are generally 
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few of the points in figure 6. For a derivation of the COBE bounds, see, for example, 
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Institute), edited by J. Chan and L. DePorcel; |astro-ph / 9 703 1 96| . 
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[37] We point out that for the run in figure 8 the production is so low that the x variance 
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negligible in comparison to the variance from resonant production. There are, however, 
a few runs with very little production, and for these the variance we record on the 
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2.4 X lO^^GeV. 
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FIG. 1. (a) : The stability chart of the Mathieu equation, Eq. (|2|). The dark regions correspond 
to stable solutions while the light regions correspond to exponential instabilities. We also show 
contours of constant fi, where /x is the instablity index. The curve fi = divides the parameter 
space into stable and unstable regions. The contours shown are n = {0, 0.1, 0.2, 0.3, 0.5, 1.0, 2.0, 
3.0}. The lines Aq = 2\qQ\ and Aq = — 2|go| are also plotted. (The plot was generated numerically, 
and in order to keep the file size small we used a fairly coarse grid. The instability bands really 
extend all the way to go = at the points = n- = 1,2,3,.... Also, for negative the 
narrow regions of stability form continuous bands rather than the "island" structure shown.) 

FIG. 1. (b) : The stability chart of the Mathieu equation for large Aq and qq. As in figure 1(a) 
the dark regions correspond to stable solutions while the light regions correspond to exponential 
instabilities. The contours shown are = {0, 0.1, 0.2, 0.3, 0.5, 1.0 }, and the line Aq = 2|go| is 
also plotted. Notice that the distance between resonance bands above Aq = 2\qQ\ for fixed qq is 
6Aq ~ I'i'ol^''^- Notice also that the instability index ^ decreases rapidly with increasing Aq. 

FIG. 2. The expectation values of the fields as a function of time for negative coupling, |go| ~ 35, 
r = 10 (m^ = = 0, = 10"^^ = 10"^ g = -10"^°). 

FIG. 3. The expectation values of the fields as a function of time for positive coupling, qq 35, 
r = 10 (m^ = m^ = 0, = IQ-^'^, = 10"^ g = 10"!°). 

FIG. 4. The variances of the fields as a function of time for negative coupling, |go| ~ 35, r = 10 
(m^ = = 0, = 10-12, A^, = 10"^ g = -lO'iO). 

FIG. 5. The variances of the fields as a function of time for positive coupling, qq ~ 35, r = 10 
(m^ = = 0, A^ = 10-12, A^, = 10"^, g = IQ-^^y 

FIG. 6. (a): Peak and valley lattice variances at the beginning of the scattering regime as a 
function of |go| for negative coupling and two values of r. 

FIG. 6. (b): Maximum physical peak and valley variances reached in the expanding universe 
as a function of |go| for negative coupling and two values of r. 

FIG. 7. (a): Maximum variances as a function of for several values of {qqI, with r = 10 and 
g < 0. The inflaton is massless (m^ = 0) and A,^ = W~^'^ throughout. 



FIG. 7. (b): A blow-up of the qq = 



—350 curve shown in figure 7(a). 



FIG. 8. (a): The variances for a set of parameters for which there is httle particle production. 
\qo\ = 350, r = 10 (m^ = 0, = 2.9 x lO'^^M^, = IQ'^"^, = 10'^, g = -10"^). 

FIG. 8. (b): The x field occupation numbers at various times for the run of figure 8(a). 

FIG. 9. (a): The variances for a set of parameters for which there is significant particle pro- 
duction. \qo\ = 350, r = 10 (m^ = 0, = 4 X lO-^^M^, = 10"^^ = 10-^ g = -10"^). 

FIG. 9. (b): The x field occupation numbers at various times for the run of figure 9(a). 

FIG. 10. Maximum variances as a function of for |(j'o| ~ 350, r = 10 and g < 0. The inflaton 
is massive (m| = 7.2 x lO-^^M^, A^ = lO'^^, A;^ = 10-^ g = -10"^). 

FIG. 11. (a): The x field variances for three runs with slightly different x masses 
{ml = 7.2 X 10-13m|, A^ = 10^^^ A^ = 10-^ g = -10^^). 

FIG. 11. (b): The evolution of the x field zeromode xo for the three runs with slightly different 
X masses shown in figure 11(a) (m^ = 7.2 x 10-^^M|, A^ = 10"^^ A;^ = 10"^, g = -IQ-'^). 

FIG. 12. (a): Maximum peak and valley variances as a function of \g\ for = lO^^GeV. The 
inflaton is massless (m,/, = 0), A,/, = 3 x 10~^^, g < and A^^ is adjusted to keep r = 10. 

FIG. 12. (b): Maximum peak and valley variances as a function of 1^1 for = lO^^GeV, 
niff, = lO^^GeV, and A<^ = 3 x 10~^^. Here A,^ is adjusted to keep r = 10 and g < 0. 
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